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5.3 ESTIMATION QF REGRESSION CO-EFFICIENTS 
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estimated at the 42 distances by taking appropriate average 
of the products of the object function values, thus obtaining 
70 values for 42x1 vector of covariances. 

Using the simulated projection data for 100 views and 64 
rays, reconstruction of these 70 images were done by CEP 
algorithm with the ’Sine’ window given by 
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Tables 5.4 - 5.11 give these errors £or samples 2-9 

respectively, while Tables 5.12- 5.19 give the normalised 

errors in the same order. These tables show that ESE is 

reduced when either R is incresed or when the covariance 

c 

kernel becomes smoother. The rate of decrease in error with 
is faster for images having smoother covariance kernels 
than for those having comparcLtively less smooth covariance 
kernels. Fig. 5.1 and 5.2 shows the actual simulated images of 
different samples and Fig. 5. 3 to Fig. 5. 8 shows the 

corresponding reconstruction with R^ equal to 16,32, and 64. 
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CHAPTER 1 


INTRODUCTION, NOTATION AND LITERATURE 


1.1 INTRODUCTION AND OUTLINE 


Computerized Tomography (CT) is a technique of 
reconstructing a function (representing the cross-section of 
an object) from its line or plane integrals called 
projections. This technique has been successfully applied in 
the fields of medical imaging, nuclear engineering, radio 
astronomy, electron microscopy. X-ray crystallography and 
geophysics etc. Continuing research for developing more 
efficient algorithms to improve the quality of reconstructed 
images calls for a theoretical analysis of various errors in 
the process. 

A large number of algorithms for CT reconstructions have 
been developed, one of them being the Convolution 
Backpro j ection (CBP). The CBP algorithm is probably the most 
widely used of all methods, as it la fast, efficient, 
reasonably accurate, and easy to implement on computer. 

Non-statist ical aspects of the analysis of error in CBP 
have been extensively studied by Natterer (1980,1986) and 
Munahi ( 1988) . This dissertation continues this study further 
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TABLE 5.11 (Contd ) 


PIC 

R 

c 

ERl 

ERS 

ERMax 

rw 

E 


16 

0.90792 

2.10333 

9.50445 

2 . 18561 

6 

32 

0.38250 

0.26267 

5.31303 

0.39268 


64 

0.22571 

0.06583 

1.38080 

0.07544 


16 

0.91043 

2.06865 

9.26359 

2.16035 

7 

32 

0.37887 

0.25974 

4.90167 

0.39173 


64 

0.21988 

0.06311 

1.36921 

0.06629 


16 

0.90836 

2.09990 

9.12071 

2.18624 

B 

32 

0.38016 

0.25937 

4.75157 

0.40057 


64 

0.21837 

0.06216 

1.34444 

0.07165 
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by developing some statistical error estimates for the CBP 
method . 

Ifi order to make the presentation as self contained as 
possible a review of the relevent literature is done in this 
chapter (section 1.4). In section 1.2 we introduce 
preliminary definitions and notation. The CBP algorithm is 
explained in section 1.3. 

In chapter 2, the problem is formally stated and 
various assumptions made to develop the theory are given. 
Various classes of spatial stochastic processes are defined, 
which play important role in the analysis of expected squared 
error (ESE). 

In chapter 3, a modulus of continuity of the covariance 
kernel is defined and an estimate of ESE is given in its 
terms. An asymptotic formula giving the order of the ESE as a 
function of cut-off frequency is also developed for those 
spatial processes, for which the covariance kernel is very 
smooth in the neighbourhood of;^origin. In section 3.5 we give 
an estimate of the ESE for the class of spatial processes 
which have very smooth covariance kernels. 

The main theorem of the thesis is proved in chapter 4, 
where it is shown that the spatial processes (cross-sections) 
can be characterized by the order of the expected squared 
error in reconstruction by CBP when the convolving functions 
belong to a certain class. 
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16 
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64 

0.05889 

0.02481 

0.01464 

0.09315 

0.03743 

0.01874 

0.14176 

0.02547 

0.00489 
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0.05921 

0.09266 

0.14050 
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32 

0.02464 

0.03733 

0.02548 


64 

0.01430 

0.01840 

0.00431 
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Chapter 5, suiBEiarizea the nuiaerical results of 
reconstructions performed using a discrete implementation of 
the CBP algorithm. A linear relationship has been established 
between the error in reconstruction by CBP and the covariance 
kernel of the process, where the covariance kernel is 
estimated by sample observations (images). The teat cases 
have been simulated as random images from certain classes of 
spatial stochastic processes on computer. 

Chapter 6, deals with a related problem of obtaining 
statistically optimal convolving functions in discrete 
implementation of CBP. A method for derivation of a 
convolving function is given which optimally takes care of 
(i) the process giving rise to the object image to be 
reconstructed, (ii) the data collection geometry, (iii) the 
cross-sectional discretization scheme, and (iv) the 
projection data noise, while at the same time it minimizes 
the ESE in reconstruction. 

1.2 DEFINITIONS AND NOTATION 


Ue will denote by k", the n-dimensional Euclidean space. 
Ue will restrict ourselves here to the case n=2 . denotes 

the positive half of the real line JR, i.e.,D?^=^x -idR :x ^ Oj 

„ \t.Srt'€S,eV c'5 „ 

The set of directions in K axie t4ie unit vectors in ER ^ 

they form the set S””*= x « CR^ : ||x|| = 1^, where [jxjj is the 
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I£ Q. is riofi-sifigulat- , it follows that 




and 

thus 

^n 

provides a 

strongly 

consistent 

estimator 

of 

■ 

If a 

is 

not known 

to be non-singular, 

any of 

the 

following 

two methods may 

be used to 

estimate q 

from 

this 


(I) TRUNCATED S.V.D. APPLICATIQM PROCEDURE : Let the rank of 

matrix 0 be r ( < K), Consider the singular value 

decomposition (S. V. D, ) 
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where > . . . S: are singular values of and , 

A, It 

and V are the unitary matrices 
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Euclidean length of the vector x « 

Let = -^ ZCu) ; u « K* J be a spatial stochastic 
process. In our case a realization z(u) of ZCu), u K*, 

corresponds to an observed image (cross-section). Thus, ^ 
denotes an ensemble of all images belonging to a class, z is 
a real valued function on iR* representing the sou^-after 
density distribution or object function. In section 1.3 and 

1.4 we shall use the notation f in place of z to denote the 

object function. Ue will restrict z to take nonzero values 
only inside a unit ball in 0?’^ (or D?*) centred at origin and 
the support of z (or f) will be denoted by fJ. 

Parallel Beam Geometry [Fig. 1.1] : A line in can be 

specified by two parameters; its (signed) distance s from the 
origin and the angle & it makes with the y axis. Ue denote by 

p(s,6) that function of two variables whose value for any 

\V-t 

(s,©) is defined as>line integral of z along line (s,6>), i.e. 


PCs, (9) 


z( s cosS -t sin6, s sin& + t cos^) dt 


where the limits of integration depend, in general, on s,& 
and O. For any s and &, (s,S) and (-s,e+fi ) represent the 
same line in 0? , so that p(s,6») = p(-s, S+n). 

The collection of values of p(s,^) for various s and & 
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TABLE 6.4,: ERROR IN PETCATSTRIP IKTEGRAL(NOISE-FREE) CASE 



IQIH 


RKS 

SHE 

CSI 

LIN 

3.715 

4.006 

3.909 

4.064 

4.277 

FAC 

3.630 

4.455 

4.325 

4-499 

4.760 

BRA 

4.964 

5.052 

4-946 

5.088 

5-248 

STA 

3.446 

3.590 

3.506 

3.626 

3.759 

COR 

14.886 

14.955 

14.909 

15-000 

15.046 

SAT 

4.295 

5.401 

5.200 

5.501 

5.828 

THO 

4.272 

4.564 

4.451 

4.629 

4-857 

SMO 

0.955 

1.026 

1.026 

1.026 

1.055 

PLA 

3.141 

3.282 

3.229 

3.299 

3.440 

BOY 

3.926 

3.867 

3.750 

3.926 

4.140 

MON 

2.761 

3.265 

3.177 

3.330 

3.506 

UOR 

3.554 

3.692 

3.649 

3.713 

3.771 

CHA 

6.110 

6.964 

6.777 

7.044 

7-310 

FIN 

3.105 

3.670 

3.589 

3.710 

3.811 

CRO 

3.824 

4.014 

3.928 

4.048 

4-169 

JET 

4.385 

4.535 

4.442 

4.573 

4.741 

GEO 

2.812 

3.203 

3.105 

3.276 

3-472 


AVERAGE 4-340 4.679 4.583 4.72? 4.894 
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ia known aa projection data and thia method of data 
collection ia aaid to be parallel beam geometry. Clearly, in 
the parallel beam geometry for data collection we require 
many aource-det actor poaitiona, which can be either time 
consuming or costly. 

Fan-Beam Geometry [Fig. 1.2]: In this case a single 
rotating source ia sufficient for an outer ring of stationary 
detectors. The projection data for this case will be denoted 
by £(£>■, /?) , where o is the angle between a data ray and the 
central ray and is the source direction. 

TbeFourier transform (FT) of a function f in K will be 
denoted by f. Thus, 


CD CD 


f CX.Y) 


J^-i2n(Xx+Yy) 

-OD 


fCx,y) dx dy 


in Cartesian co-ordinates and 


TT 00 


f (RcosS, RsinS) = 


e 


- i2nRccoa 


r j dr d^ 


0 -00 


in polar co-ordinates. The inverse Fourier transform of a 
function g will be denoted by g. Thus, 







r strip integral 
(b) RAM filter, 

, (e) original 




h lO 


7 



Fig. 1.2 Fan-beam data collection geometry. 




’0 c on s txriic tion of . Cluro^^rnosorn^s in Pii^T s'cirii^.^i' i,nt0^ir0j 

noise-free case Ca/ CSI nAJ. filler, 

(c)RKS filter, (d) Shepp filter, (e) original 
picture, (f) optimal filter. 
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g (r,^) 


fi m 

r r 

1 ' /n ^ i2nRrcos (©-<(&) ,„, 

g(RcosSf, RsinS) 6 <. T'-' j£j 

0 -ao 


The notation PCx) is used for the probability density of 
an n-diinensional random vector at x « k", and P(x) will 

denote the probability density of a random variable at x « E. 
A random vector X is said to follow a multivariate normal 
distribution if, 

P(x) = (1/ expC-| Cx-y)' S'^Cx-y)) 

where fj is called the mean vector and E the covariance 
matrix. A random variable X follows Poisson distribution with 
parameter X if, 

-X . X 
e X 

p( X = X) = — — j 

PCx/y ) will denote the conditinal probability density of 
random vector X given X ~ thus 

P(x/y) = PCX, y) / PCr)- 

1.3 THE CBP ALGORITHM 

The CBP algorithm is derived from the well known 
Projection slice theorem [ Bracewell (1956) ]. Let p( R, S) 
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denote l-dimensional FT o£ pCs,^) with respect to the first 
variable s, and fCRcos^, Rsind) denote the FT of f(c,^^)in two 
diiaensiona, i.e.. 


pCR,^) 


Or 

p(s,d) 


-i 2nRs 

6 


ds , 


-00 


and 


n a> 


f (RcosS, RsinS) = 


js -i2fiRrcosCS-(^) I , j .. 

fCr, 4 ^) e ^ |r} dr d^^ 


0 -00 


Then this theoreai states, that 


p (R, &) = £ (Rcos©, Rsin©) - ....(1.3.1) 

A two dimensional inverse FT of (1.3.1) leads to the well 
known tomographic formula. 


n 00 


£ (r.-jfr) = 


p (R.e) ^^2nRrcos(& 


. .(1.3.2) 


0 -cot 


In practical implementation of equation (1.3.2), the factor 
jRj in the integral is replaced by |Rj U(R), where U(R) is a 
suitably chosen function called a window function. Normally, 
UCR) is an even function of R, the Fourier frequency and has 




-«v> *■* 


truction of Lincoln i: 
case (a) CSI filter, 

, (d) Shepp filter, (' 

tiraal filter . 




Mi 





hepp 

liter 
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a compact support [-R , R ]. If the frequency content of p is 

e •& 

also supported in [-R , R ], i.e., if p (R, vanishes for 

a 

iRl > R , then the filtered version with UCR) = 1, for 

* * CS 

R <s [-R ,R ], and zero for iRl > R , 

e “S * * -e 

n m 


£ Cr, 4i-) 


r ^ 

J p (R, e) 


^i2nRrcos(S-4i^) 


Rj U(R) dR dd 


0 -CO 


...(1.3.3) 


agrees with fCr, ^) as given in (1.3.2). In the general case 
introduction of U(R) with a finite cut-off R helps in 

Ci 

stabilizing the singular nature of the integral in (1.3.2). 

An implementation of (1.3.3) in^spatial domain is as 
follows: 


fT as 


f (r, ^J.) 


p (s, &") q(s' - s) ds dd 


0 -oa 


(1.3.4) 


with s' = r co8(,&-4>') , and 


q(s) 


as 

|Rj U(R) dR 


-05 


(1.3.5) 


The. inner integral in equation (1.3.4) is a one 
dimensional convolution while the outer integral corresponds 
to the backpro jection operation. Thus it is named convolution 
backpro jection method. The CBP is also known as filtered 
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backpro jection method indicating filtering of projection data 
p, fay window function U(R). The practical implementation of 
CBP as given by equation (1.3.3) was initially suggested by 
Bracewell and Riddle (1967). Ramachandran and 
Lakshminarayanan (1971) described the method which took care 
of the discrete nature of projection data. Both Bracewell and 
Riddle and Ramachandran and Lakshminarayanan used the band 
limited window 


1 


U(R) = 


R < R 

e 


(1.3.6) 


to jRl > R^ 

in their work. Shepp and Logan (1974) suggested the use of 
’sine’ window given by 


W(R) = 


sin C ^ 


c nR /2R^ 5 


R < R 


0 , jR| > R^ 

Use of several other window functions has been reported in 
V^.’e.l it erature . Ue note that these window functions are radially 
symmetric in the Fourier domain. 

In practical situations where projection data is 
available for only a finite number of rays and views, a 
recommended value of R , based on sampling theorem [see Jerry 
(1977)] is given by R ^ l/2As, where As is the spacing 
between two^ consecutive data rays. 
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1.4, A REVIEW OF THE LITERATURE 


Ift this section an attempt is made to give a reasonable 
review o£ the available literature in reconstruction methods 
and their error analysis. A reasonably complete bibliography 
of introductory literature in CT can be found in Deans (19S3) 
and Herman (1980). These books also give an excellent survey 
on the applications of CT. A state-of-the-art review is 
contained in Harch 1983 issue of the Proceedings of the IEEE, 
which is devoted to CT. A good account of the mathematics of 
CT is available in Herman and Natterer (1981), Natterer 
(1986), and Stark (1987). 

The first mathematical formulation of the problem of CT 
is due to Radon (1917). He defined an integral tranform R, 
which maps a function f on k" into the set of integrals over 
the hyperplanes of k". Kore specif ically^if & , and s e 
0? , then 


(Rf) (©, s) 


f(x) dx 


f(s & + y) dy 


X.& = 3 


d 


X 


where = •£ x « : x.S =0 }■, and x.& denotes the inner 

product of the vectors x and §. This transform is known as 
V^e. Radon transform. Louis and Natterer (1983) defined the x-ray 
transform P, which maps a function on k" into a set of its 
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transform, SIAM on hhm^x-icctl 

pp. 456-467. 
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line integrals, i.e., if & <e and x e , then 

<u 

CPf) Cf, X) = fCx + t 0) dt. 

-a, 

For n = 2, with x S'*", p and R coincide, i.e., 

CR£) (S, s) = (Pf) (s.e, X). 

Another way to label the straight lines in D?*' leads to the 
divergent beaai tranaforiaation 

OCs 

(Df) (x, S) = fCx + tS) dt 
0 

where x is thought of as the source of a ray with direction 
&. In emission tomograpy the transform which decribes the 
data is known as the attenuated Radon transform 

CR^f) (s, S) = fC X) ^ dx . 

X .S=s 

where fj is the attenuation coefficient and f is the activity 
distribution . 

All these transforms describe the generation of 
projection data using different data collection methods. The 
data in Positron emission tomography (PET) corresponds to the 
values of the product 
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-R 

e ^ Rf 

and do^not give rise to new matheinat ical probleiaa. 

Thus the reconstruction problem o£ CT simply calls for 
the inversion of Radon transform in The first inversion 

formula was given by Radon (1917) himself. It is given by 


fCx) 



00 

dF (s) / s 


0 


where , 


F Cs) 

K 


1 


pCf, x.f + s) d§ , 


and p(s, ©) is the projection data. 

Bracewell (1956) gave the Projection theorem, which 
expresses the relationship between an image, its projection 
data, and its two dimensional Fourier transform in a way that 
is both mathematically precise and intuitively clear. In 
mathematical notation it says that 


p (R, = f (RcosS, RsinS) 


(1.4.1). 


Later, Bracewell and Riddle (1967) obtained the famous 
tomographic formula 


fCr, ^) 


TC 00 

r r 

I I p(s, ©) q(s’ -s ) ds dd 


0 -00 
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by taking a two-diiaenaional inverse FT. This is known as 
Convolution Backpro jection (CBP) method or filtered 
backpro jection method. A practical implementation of this 
method was derived by Shepp and Logan (1974). 

Herman and Maparstek (1977) derived the reconstruction 
formula for "fan-beam” geometry ( or divergent beam 
transform ). This reconstruction me 
equations : 

2n S 

(1/ sin(e>‘ '-O') ) D gCo-, /3) do d/3 

V 

0 -s 

where <5 is the object radius. 


f(r, = (l/4n“) 


thod is given by^fol lowing 


o' = tan 


r cos (/3 - <^) 

D + r sin(/3 - (ji») 


g(o, ro =-^ 


f ^ _ fg 1 

{_ ^o J 


and 


U = [ C r cos(/9 -0) > + C D + r sin(/3 - ^) } ] 


2 


where g(o, /3) is the projection data collected in fan beam 
geometry and D is the distance of source from the object 
centre. It is also shown that the singular integral in 
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(1.4.2) can be regularized in a number of waya each of which 
givea an efficient numerical evaluation of the function 
atrueture and apeed is same as that of "convolution 
reconstruction technique”. The reconstruction formula for 
divergent beam data as in (1.4.2) is also derived 
independently by Lakshminarayanan (1975) by a method 
different from that of Herman and Naparstek. 

Instead of inverting the FT in projection theorem 
exactly the use of numerical inversion with the fast Fourier 
transform (FFT) is the so-called Fourier method of 
reconstruction. This is also a transform method. Fourier 
reconstruction in CE has been suggested by Bracewell (1979) 
in radio astronomy and by Crowther et.al. (1970) in electron 
microscopy. The method is to use direct Fourier inversion in 
projection theorem. The problem with this arises in 
discretization hence in approximation of f on a regular grid 
from the known values of f on a polar grid. The critical 
point here lies in interpolating f in s directions. Bracewell 
(1979) has observed that failure of the standard Fourier 
reconstruction method results from interpolation. It is also 
shown that interpolation which is correct in the sense of 
sampling theorem requires a huge number of operations. Low 
and Natterer (19S1) give some optimality criteria for 
interpolation points, [ see Stark et.al. (1981a) ]. Stark 
et.al. (1981b) replaced the nearest neighbour interpolation 



17 


by a fiiot'e aophiatieated Interpolation named near exact 
interpolation method and ahowed that quality of image is same 
aa that obtained by CBP method, but time taken is less than 
that by CBP. Natterer (1985) also discusses errors in 
interpolation. 

A different type of inversion formula is produced by 
series expansion method in terms of special functions. This 
approach differs fundamentally from transform methods aa this 
is based on discretization of the image domain prior to any 
mathematical analysis, whereas in transform method first the 
formula for inversion of Radon transform is obtained and then 
dicretisat ion is done. These methods are also known as 
algebraic methods, iterative algorithms, or optimization 
theory techniques. Censor (1983) has formulated a fundamental 
model in this approach to image reconstruction problem as 
follows 

A Cartesian grid of square picture elements, called 
pixels, is introduced into the region of interest so that it 
covers the whole picture that is to be reconstructed. The 
pixels are numbered in some agreed manner, say from 1 to n. 
The X-ray attenuation function is assumed to take a constant 


value X 

i 

througliout the j-th 

pixel , 

for j=l , . . . 

,n. 

Source 

and 

detector 

are 

aastimed to be points and rays between 

them 

are 

assumed 

to be 

lines. Let a.. 

'-i 

denote 

the length 

of 

i-th 

ray 

segment 

that 

is intersected 

by the 

j-th pixel. 

for 

i=l, . . 

- 
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j = Thua phyaical meaaureiaent o£ total attenuation o£ 

the i-th ray, denoted by g. , represente the line integral of 

1 . 

the unknown attenuation £unction. Let R.£ be the line 

1 

integral of fCr,^) along i-th ray. Then g^ is an 
approximation o£ R^£. The whole model is described by a 
system of linear equations, 

n 

2 x. a. . = g. i = l m, 

t »- 4 1. 

j=i 


Algebraic reconstruction technique (ART) was first published 
as a reconstruction algorithm in a paper by Gordon et.al. 
(1970) and was later recognised to be identical with 

Kaczmarz ’ s(1937) algorithm for solving system of linear 
equations by Guenther et.al. (1974). Realizing the difficulty 
in handling such a large system of equations Gordon et.al. 
(1970) suggested an iterative process, starting with an 
initial approximation x*^ « E*' to the image vector. In an 

ic 

iterative step, the current iterate x is corrected to give 
new iterate x by taking into account only a single ray, 
say the i-th and changing only the image values of the pixels 

which intersect this ray. The discrepancy between the 

fl k 

measurement g. and pseudo-projection data >■ ^ a.. x. 

£ij-l 3 

k 

obtained from current iterate x is redistributed among the 
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pixels along the whole o£ i-th ray in proportion to a. .’s. 
Thus an iterative procedure ia given which in practice stops, 
but unf'ortunately not iniueh is known about how "good” the 
iterates themselves are in comparison to the limit. Also in 
practice limits are not achieved and one has to make 
approximations. Hounsfield (1973) has obtained ART 

independent of Gordon et.al. (1970). Censor et.al. (1983) and 
Eggerment et.al. (1981) have studied the behaviour of the 
algorithm for the inconsistent system when relaxation 
parameters are allowed. Tanabe (1971) found that the ART 
works well for both singular and non-singular systems and 
also gave the affine space formed by the solutions if they 
exist. A complete tutorial of ART algorithm is given in a 
paper by Gordon (1974). 

Based on statistical nature of data and image, certain 
statistical reconstruction methods are also proposed in 
literature. Hanson (1987) considers the projection data g as 
a weighted integral of the object function f with weight 
functions h. , which are charact erst ic functions of lines or 

I- 

strips in the plane. Thus 


CO CO 


e. = 


X, y) f( X, y) dx dy 


-ao -ao 


give rise to the linear system 
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g = 5 f - 

Uith noise in data the model becomes 

g = H f + n , 

wher-e n is the noise vector. Suppose that n follows a 
multivariate normal distribution, with mean zero and the 
covariance matrix so that the conditional probability of 

g given f, PCg /f ) is multivariate normal with mean Hf , and 
covariance matrix R , i.e., 

P (g /f ) « exp [ Cs - Hf >' r"^ Cg - Hf > ] ... Cl. 4. 3) 

If the a. px-tax-t probabilty distribution for the ensemble of 
images, is assumed to be a multivariate normal distribution 
with mean value f. and covariance matrix R , i.e., if 

-0 f 

PCf) « exp [ -| cf - fg)’ Rj* Cf - ]' Cl. 4. 4) 

then the best estimate for the object image f is that 
particular image f which maximizes the a. post&x-iax-i 
probability 

P Cg /f). PC i ) 

P Cl /g) = Cl-4.5) 

PC g ) 

Under these assumptions Herman et.al. C1979) have shown that 
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maximum a. pOBtsr-ia:iri (MAP) solution satisfies the MAP 
equation, 


r/ (f.- f) + H’ R"‘ (g - Hf) = 0 

I — y — — y.% ~ 


... (1.4.6) 


Herman and Lent (1976) have given a computer implemenation of 
Bayesian (or MAP) approach as an iterative procedure, the 
iteration scheme was given with initial guess f^, i.e.. 


= Iq- R“‘ (g - Hf*"), 

k k, k, k, k k . k 

c =£s /s s, s =(I+R^H R^H)r,J 


(1.4.7) 


where the vector r is the residual of the equation (1.4.6) 

k 

(multiplied by R^ ), the scalar c is chosen to minimize the 

k'4*^ 

norm of r , and I is the identity matrix. Similar scheme 
for Image restoration is proposed by Hunt (1977). Trussel and 
Hunt (1979) redefined the iterative scheme given by Hunt 
(1977) to improve the rate of convrgence. Hanson (1987) found 
that the scheme (1.4.7) does converge but takes at least 10 
to 20 iterations. When there is an absence of a pr-ior-i 
information, as is often the case, then a. px-icfx-i probability 
distribution P(f) does not play any role in (1.4.5). The MAP 
approach then reduces to that of maximizing the conditional 
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probability (1.4.5), or equivalently to that o£ minimizing 
the so called chi-square 


/ = ( g- Hf) 'R^‘c g - Hf). 


and then the method is called the minimum chi-aquare method 
or the least squares method. The solution of this is also 
refered to as maximum likelihood solution as this term is the 
likelihood function for normal distribution. Rockmore and 
Mackovaski (1977) have used maximum-likelihood (ML) method 
for reconstruction without the normality assumption. The 
method is summarized as follows: 

Let f(x,y) represent the linear attenuation co-efficient, and 
the incident process (as a source of photons is incident on 
object ) has rate then the rate of detected process is 
given by 


X = exp C 


j f(x, y) dx dy ) 


where is the line by which photons are detected (Collimator 
is alligned with the line JS ). It is assumed that the 
detectors observe a set of random counts, each having 
a Poisson distribution, corresponding to a series of N 
projections over time interval of length t. Let these 
projections be indexed in an arbitrary way. By virtue of the 
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properties o£ the Poisson process, they are independent, thus 
the observables, denoted by ot>ey the joint Poisson 
probabilty density function. 


N 

P Cnj^= , n^= k2 , . . . . > = FI ^ ^ 

i=l 


1, ^ , k. 

6 J. C > t 

= __ 

L 

where the are given by 


= Xq exp C 


1 


f (x , y ) dx dy ) 


1 

being the i-th source-detector line. In order to apply ttes 
the ML estimation technique the function f (x, y) is expanded 
on an arbitrary complete set of basis functions { ^^^Cx, y) }, 
and for a finite number L, it is given as 


L 

f (x, y) = 5 (X, y) = f ' ^ Cx, y) ..(1.4.8) 

I =1 


where, 

f = f ) . and ^ (x, y) = C (x, y) , . . . ,4'^ (x, y) ) 

» 2 la **** i Mmd 


Then the parameters to be estimated is the vector £. Using 
(1.4.S) it is shown that 
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r 

1 . 

L 


f (x, y) dx dy = 2 f^^ 



i =1 





(x, y) dx dy 



¥ 


i i 




Uith these notations the pcoblein is to find the beat estimate 
of f (ML sense ) given a particular set of realizations 

(k^.k^ k^) corresponding to number of counts for each 

projection where the number of counts are independent Poisson 
variates. Thus the ML estimate is that f which maximizes the 
conditional probability density function P (n = k /f ) , i.e, 
the estimate is the vector f solving the equation 


S 


P c 


n 


= k j 


i > 


0 


for transmission case the ML equation is 


N 


2^ 





V, 


U 


k. 3 

■L it 


0 

I = 1,2,..., L. 


The solution is given by 

I = C + y, 

where Q = C -In C (k^^ + z^ ) /X^t ) -InC (k^+z^ ) /X^t ) ), z is 

any vector in null space of and y a vector in null space 
of ^ Rockmore and Macovaski (1976) used ML method for 
emission tomography using the rates 
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f (x, y) dx dy 




instead of exp C 


fCx, y) dx dy y 


£ 


the approach remaining the same. The ML equations now become 
N 

2 I - V',, t + k. /£'¥',)] =0. 

L = 1,2, . . . , L. 

and the solution is given by f = (1/t) C w k, where, 

If = C >. 

Shepp and Vardi (1982) have used maximum-likelihood 
estimation in ease of Positron emission tomography (PET). They 
give a general mathematical model for PET where an unknown 
emission density X is to be reconstructed from the number of 

t 

counts, n (d) in detector unit d, d = 1,...,D. The unobserved 
emissions n(b) in b-th pixel are assumed to be independent 

Poisson variables with unknown means X(b), for b= 1 B. 

Suppose each emission in pixel b is detected in detector unit 
d with probability p(b, d), d= 1,...,D with p(b, d) a one 
step transition matrix, assumed to be known. For estimation 
of X(b) the Expection Maximization (EM) algorithm of Dempster 
et . al . (1977) is used. The algorithm contains two steps the 
expection step and the maximization of likelihood step. It 
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starts with an initial estimate X and gives a new estimate 

estimate iteratively as follows: 


, Ffft V . , . 

X Cb) 


''•old 

Cb)J. 

d=l 


n Cb) pCb, d) 

9 

J x^^^(b') PCb' ,d) 
b'=l 


b=l B. 


It is shown that this iterative procedure converges. Vardi 
et.al.C1985) have also investigated further estimation 
methods in PET such as least squares method, method of 
moments, etc. They have compared EM algorithm with other 
methods, and observed that it gives better results than least 
squares method. Lewitt and Muehllehner (1986) give 
accelerated version of EM algonbithm which takes about half 
number of iterations than that taken by of standard EM 
algorithm. 

Another statistical method of reconstructing an image 
from its projection data is known as a Maximum Entropy (MEKT) 
algorithm, which also uses cl pricrr-i information about object 
function. Let f^, €^, . . . denote pixel values of the image 
to be reconstructed and p. = r. / T r. , where r. is the 

X. 1. L 

average rate of emission of photons from i-th pixel. The 
entropy is defined as 


H CP3_. . - • .Pjn > 


n 


1 Pi ip* Pi 

i = l 


h /Ih 

i = l 
i = 1 


n 





The entropy depends only on the distribution o£ gray levels 


and not on the total intensity ^ f ^ . The observed data is 

i 

given by 


n 

g . - ^ A., f . + e . 

J Z Jt- V. j 

i =1 


j = 1 , . . . ,m. 


where e^ is independent random noise with mean aero and 
variance & .Let 


n 


m 


Q C . 


> 


i 1 « I \ - 

j=l i=l 

T S - e !d 




2 2 : 

where D =diag [ l/<y ] . With these notations the 

1 m 

problem o£ image reconstruction becomes that o£ maximizing 
the entropy 

“ C ) 

subject to constraint 

Q C f 5 = i«/2 and T f. = t Ct>0). 

to 


Zhuang et.al.(19S7) give a differential equation approach 
solve this problem. The differential equations are 


as 
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where , 

J (f: X ) = H C£^ £^) 

and ^ 

7 J = - [^ logfj,..., logf° j % C ^-1) h - X a 'D CAf° - g) 
a discretization of these equations is given by 


= exp (,p - 1) h, 

k k+ 1 

k 3: 0: C F + X, A 'D A ) f 
k - 

= C2Xj^ - X^^^) A 'D A + h + Xj^ ) A 'D g 

where = diag [ 1/f^, . . . ,1/f^ ] , X^ = 0 , and X^^^^ > X^ , and 
h is a nxl vector containing I’s. Using Gauss-;;^eidel 

technique following iterative acheme is obtained 


k k-i-1 

k > 0: CF X^f ) = (2 Xj^ - Xj^^^ 


)d + h + c Xj^) g. 


k+ 1 k 

k> 1: g =C2 -X^^^/X^> d 

+ ( 1 /Xj^) C ..., C£^ - £^'"^)/f^ y 


* f -1 > e- 


where 


L = A 'D A, 
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E = A 'D g, 

k = 0,1,2,.... 

The initial values f^, d^, d^ are as follows: 

= exp (v -1) h, 

n 

d° = Lf" = exp ( ^ - 1) C 2 Lj. y. 

i =1 



Uernecke and D’Addario (1977) give a solution of this problem 
in radio astronomy using Fourier synthesis. Minerbo (1979) 
defined the entropy as 


H (f) = 


f (x, y) log C fCx, y) A ) dx dy 


D 


where A is the area of D.The constraints are 
a 




+ 1 ^ 




f (s COS 6 f. - t sin s sin B, + t cos BJ dt ds 
4 j a 4 


a , ~c» 

jiR 


where s < 

ji 


< s. 


iMip 


HI — J •" 1, 

is a set of abscissas for j-th 


veiw- The iterative algorithm is given as 

= 1 , m « Jl , (Ji is set of m values for' which g. 2:0) 

jm - aw 
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H 


k+1 

■ jffi 





A g 


jra 


n 2 H. 




at. Cs 

hrt 


cos 

J'- 


t sin ) 

j*- 


i^j 


m e ^ , j = k rood C J + 1 ) 


ms v4l, j k mod (J + l) 


where , 


H . = exp ( X . - 1 ) / J , 

jm jm 



s . < s < s 

jm 

otherwise 


jm+1 


Wood and More (1981) give a minimtiffl variance estimator in 
CT , they assumed the measurement vector g = C fi^| ••• | 5 

and projection matrix H ' = [ H^' j ... | ] where g^ is 

the nxl measurement vector for j-th rotational position thus 
the usual measurement equation is g. = H.f + n, the minimum 

-J j- 

variance estimator for vector f is given by 

where £ and n are assumed to be random variables, the a 
priox^i information is given by 

= E (f), Pq =e cc f - Iq) cf - igy) 

and the assumptions about the measurement noise are E (n) =0, 
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R = E (n n E (n f ') = 0. R = C H H ' + R ) can be 
identified as the measurement covariance matrix. A fast 
implementation of this method on computer is also given. 

The problem of image reconstruction when projection data 
is available for limited number of angles, is known as 
"limited data problem”. It has been considered by many 
researchers. Hanson (19S7) and Medoff (1987) give a review of 
this problem in context of statistical and deterministic case 
respectively. Medoff et.al. (1983) have discussed the 
iterative convolution backproj action algorithm for limited 
data situation to reduce artifacts. Invouye (1979) gives a 
solution to this problem by Fourier reconstruction method 
using the analyticity of the Fourier transforms. Assuming 
that data are given at following L points. 

&= B , & , , d = & (0 < 0 < 0^), the expressions 

4 a n ft ^ a ft' ^ 


L M 

_ _ _ 

(a.) y I > [ a (<uj) cos2m 0. + b (6j) sin2m 0. ] - Ref (ij , 0. )j 

£j ^ j£j m ^ m i. l * 

i=l m=l 

L N 

Cb) T I y t c (i«>) cos(2n + 1) 0. + d (6j) sin(2n +1)0. ] 
i=l n=0 

-Imf (^, 0p }* 

are minimized with respect to constants a^, b^, c^, d^, where 
M and N are finite numbers for which Fourier series expansion 
is taken. The values of a , b , c , d , 

wi m Tt n 


are obtained for 
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given s which aiinimizes (a) and (b) and making use of 
these constants the missing data is extrapolated. Gilbert 
(1972) proposed a method of reconstruction known as 
Simultaneous Iterative Reconstruction Technique (SIRT). 
Lakahaiinarayanan and Lent (1979) have shown that this method 
strongly resembles the Richardson least squares algorithm. 
They also gave generalized SIRT (GSIRT) algorithm, which 
convergefiow faster than the SIRT. 

The mathematical analysis of the error for CBP method 

has been done firstly by Natterer (1980,1986). He has given 

the error for CBP algorithm for object functions belonging to 

certain Sobolev spaces. One of the main theorems of Natterer 

(1986) states that for f -e C™ (fJ) ; i.e., f being infinitely 

differentiable with comapct support ft, the discretization 

errors in CBP are given by e^^ and e^ , for some ^ , 

0<.#<1,R < and R 1/2 As, 

c c 


n ^ I C ''c 


*2 - • ^ 'l (a) ^ '■^•‘'5 • 






= 2n sup 
S e S 


n-1 


jRj jfCR, a)j dR , 


R >R 


m is an integer such that the quadrature formula for 
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integration is exact for the function belonging to H ' , the 

2m 

class of even spherical harmonics of degree 2m, and y ) C^, m) 

is such that 0 < yj (^, R ) < C e”^‘‘'^^^c; with R > R , 

c c 

C(#) and are positive numbers. These errors are introduced 

by discretisation of comntinuous integral in the inversion 

formula. If a function is essentialy R^-band limited, then 
* 

the term is negligible and error in CBP is strictly due to 
discretization of the backpro jection integral. 

Natterer (1980) has developed reconstruction error 

, ot 

bounds for the function f belonging to Sobolev space H , « 

real, with norm given by 




(1 + I R|® )"l f (R) j ^ dR < 


m 


The r econstruct ion error r is defined as. 


rCa, n, ?) 


sup [ II f^ -f^llj^ 


(Cl) 


«^''h«(Ci) ^ 1 


II Rf, - gll < ? , 
i=l,2 for some g J 


The lower bound of 
r(a, n, ? 


r is shown to be given by 

) i =C«, f!) 1- 


where c is a constant and projection data p is known for n 
lines with a root mean square error ? , Rf^^ is the 


Radon 
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tf-ai'ia£of-ffi Of tfue pfojeetiori data for f. 

Nattsrer in thia atudy of errora haa aaaumed that the 

crosa-sect ion (image) haa band limited projection data. But 

in general object croas-aect iona have a finite aupport, hence 
the projection data is necessarily not band limited. 
Nevertheless, in the implementation of CBP method, the 

approximaton obtained is band limited. Munshi (1988) has 
studied the error in CBP algorithm due to thia band limiting 
aspect. The analysis assumes perfect (noise-free) data. 

The total error, e , in reconstruction is assumed to consist 

T 

of three additive parts: 

e = e + e + e , 

T R B e 

where, e is the inherent error due to finite cut-off R , e , 

R e B 

and e are discretization errors due to discrete 

* 

implementation of backpro ject ion and convolution integrals 
respectively. The error estimates for the inherent error due 
to finite cut-off R are obtained. The estimates are given 

■St 

for f belonging to certain circular Sobolev space O and a 

2 to 

characterization of functions f is given in terms of 
interpolation space ^ , where m is a positive integer and 

am 

0 < /? < 2m. It is assumed that the convolving functions 
satisfy certain regularity conditons. 

Bertro et.al.(19SQ) considered the inversion of Radon 
transform as ill-posed problem in estimating f from the 
available data g = Rf + h, where h is the noise in the data. 
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The error in estimation is defined as 

I £- f Ip 

where F and G, solution and data spaces, are considered as 
the spaces of square integrable functions and f has support 
on a unit circle with following constraints 

Ci) I Rf - g < ?, and Cii) I Bf |p < E 

with B : F — *■ F, a constraint on f and E is a prescribed' 
constant. The stability estimate of this error is given as 

6(^,E;g) < 6C?,E) 

The estimate is given by 

Sa , E) = y2 ^ [ |vl"^ + C? /E)^ 

Madych(1980) also has defined Sobolev spaces (D) as 

the class of those Schwartz distributions, f, on ER which are 

P 2 

supported in D with the property that f - f is in L (IR ) 

y 

for all y and £ denotes the y translate of £, with finite 

norfii 

M = in>o i 1=^1“ « 'y - ‘ Il^c k“) }■ 

Then he gave a reconstruction formula R^, using n projection 
data such that for f being in L^CD) 

\ p£ (x) = 0^* f Cx) .-.C1.4.S) 

where n is the total number of lines for which projection 

2 

data is given and ~ (!/") *!( [xj ), with 
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-2 2 

q(t) = c(t-f.) (P (t)) where P is the Legendre polynomial 

1 ni m 

o£ degree m. Then the theorem states that: Suppose £ is in 
P 

(D), 0<a:^l, l<p< 0 £>. Then given n distinct projection 
data p and the reconstruction formula (1.4.S),it follows that 


I f 




C 


n 


f l^p where c is a constant 
•a 


independent of £ and n. 

Latitach(19S9) gave an error analysis for spline 
inversion of Radon transform in the set-up of Sobolev spaces 
used by Natterer. 



a. i SOME DEFINITIONS AND ASSUMPTIONS 


In this thesis we assume that the object function ZCu) 

jg 

is a random variable for each u « K . Thus realization z of an 
object cross-section is a function on iE to K or K = [0, eo ) 
and = {ZCu): u -s 0? } is a spatial stochastic process 

representing the ensemble of all cross-sections of interest 
of which the given cross-section is a sample realization. 

DEFINITION 2.1.1 Translation of an image ( cross-section ) 

2 

ZCu) by a constant c = (c^, c^ « K is denoted by 

T^ ZCu) = Z Cu + e ) 

DEFINITION 2.1.2 Rotation of an image ZCu) by an angle & 

2 

about a point x e K is given by 




"^1 ^ 

COS 

& + 

Cuz 

- Xz ) 

sin & 

R ^ZCu) = Z 

X, & - 


- Xj^) 

sin 

e - 

Cuz 

" ^2 ^ 

cos & 


DEFINITION 2.1.3 The expected value function of the spatial 
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stochastic process ^ is given by 


E ZCu) 


Z (u) dP (Z) 




and the product Kioment function is given by 


K (u, V) = E [Z Cu). ZCv) 3 


Z (u) Z Cv) dP CZ) 




where P is a probability measure on the space ^ 


DEFINITIOM 2.1.4 

(u, V ) = E[Z (u) -E Z(u) ][Z Cv) - E ZCv) ] 
is called the covariance kernel of the process 


DEFINITION 2.1.5 The process ^ is said to be invariant under 


translations if 


and 


Z « implies that T Z e ^ 

c 


PCT^Z«B> = PCZ«B>, 


V c s and V Borel set B in £R* 


DEFINITION 2.1.6 The process ^ is said to be invariant under 


rotations if 
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2 e implies that R e ^ 

X , ^ 

and 

P<R^ *Z«B> = PCZ«B>,V«. X. ^ 

DEFINITION 2.1.7 A spatial stochastic process ^ is said to 
be second order wide sense (or weakly ) stationary if 
CH» Y ) depends only on the distance between the points u 
and V 

In practice the object function Z(u) is composed of two 


parts: m(u) = E ZCu), 

the fixed part 

and 

the random 

part 

IC^), i.e- Z(u) = 

+ S(u). The error 

analysis 

for 

the 

fixed part 

CBP algorithm 

has 

been 

done 

by 


Munshl (1986). Here, we consider the random pai't , 

S(ti) = Z(u) - m(u) only. Thus without loss of generality 

2 

we assume that E Z(u) = 0, V" u « K . Also then, the product 
moment function is the covariance kernel (u, v). 

THEOREM 2.1.1 : Let ^ be a spatial stochastic process 
invariant under translations and rotations. Then if E Z(u) 
and Z(u). Z(y) > exist, the process ^ is second order wide 
sense stationary. 

PROOF: He have to prove that (u, v ) is a function only of 
the distance d =| u - y j[ between the points u and y. Now, 
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V) = E C ZCu). ZCv) ) 


ZCu). Z(v) dP(Z) 


zfa.)- R„, aZ(v) ap (R^_ 
sc 


Zn 


Z(u) 


•( 


Zn 




ZCu). Z^^ ^ ZCv) dPCZ) 


ZCu +c). Z^ — ZCv + c) dP 


ZCu'). Z^^ - ZCv') dPCZ) , 

jar 


where u' = u +c^, v’ = v’ + c^. Now d = | 


dPCZ) 


CT Z) 
c 


U.’ ■ V’ 1. i-e. 
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distance between points u’ and v' . Hence (u, v) is a 

function of d only. 

2.2 THE ERROR IN CBP ALGORITHM 

The discrete implementation of the CBP algorithm as 
given by Levitt (1983) involves three kinds of error. . 
Firstly, the object function Z is supposed to be a sample 
realization of the process and the reconstructed image Z 
by CBP algorithm is also of statistical nature. Thus a random 
error Z - Z in reconstruction is introduced. A measure of 
this statistical error in reconstruction is the expected 
squared error (ESE) E C 2 - Z = eJ,. 

Besides the inherent statistical error there are two 
types of discretization errors as explained in the following: 

Suppose that the image is discretized on MxM pixels, and 
the projection data p(s, S) is available for K rays and L 
views. Let As and A& be the ray and the angular spacings 
respectively. Then 

L 

Z (mAx^, nAx^) ^ AS J P (s', ) , ... (Z.2.1 

1=1 

where, s' = m Ax. cos + n Ax_ sin d . The values of 

1 I 2 1. 

PCs’, the convolved projection data, are required for 

each of HxN points (mAx^, nAx^), 0<m<M, 0<n < N. 


For 
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large number o£ rays, a practical approach is to evaluate 
p (kAs, for K k and then use ’’inexpensive” 

interpolation to estimate the required MxK values of p from 
these K calculated values. Thus the convolution in (1.3.4) is 
approximated by two operations on discrete data, a discrete 
convolution for K £ k £ 

p^Ck 'As, > = (As) J p CkAs, ) q CAs (k ' - k) > 

k=K‘ (2.2.1) 

followed by an interpolation operation 

PjCs' , = (As) 2 p^Ck 'As, ICs'-k'As^ ...(2.2.2) 

k’ 

where I(s) is an interpolating function and the number of 
terms in sum depends on the width of non-zero part of I(s). 
For example, the mid-point interpolation takes only two 
points for summation. In practice popular interpolation 
schemes are the linear interpolation and nearest neighbour 
interpolation schemes. 

Ue note that the discrete convolution makes use of 
equispaced projection data. 

So the errors in reconstruction by CBP at a point x « K 
are given as:— 

(i) The expected squared error (ESE) denoted by E^ : 

(X) = E^ CZ (X) - Z (x) 
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where Z(x) is same as £ in (1.3.4) with Cartesian 
co-ordinates (x^^, x^) instead of polar co-ordinates (r, i^) . 


The 

ESE 

E® 

depends 

only on 

the 

f init e 

Fourier frequency 

cut- 

■off 

R . 
c 

1 1 can 

be made 

zero 

in the 

case of R -band 
c 


limited projection data with window function UCR) = 1 fot' 

R e [~Rg» SLiid UCR) = 0 for" |R| > R^- Our analysis in this 

thesis is for this error . 

(ii) The convolution error E^, which is due to discrete 
implementation of the convolution integral In which 
interpolation also takes place. However, this error can be 
made arbitrarily small by taking a large number of rays. 

(iii) The Backpro jection error Eg, which is due to 
discretization of the backpro j ect ion integral. This error 

can also be made arbitrarily small by taking large number of 
views and rays. 

So the total error in the reconstruction at point x « tR, 

using CBP algorithm is a combination of these three errors 

2 2 

E and of which E ^ is inherent in the process while 

C Jd ^ 

E and E_ can be kept under control, 
c B 

2.3 WINDOW FUWCTIQNS 

The window functions considered are of the form 

CR ) = W C R / A ) 

where U is a fixed function and A is the same as R^. 
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Let 


CO CO 

K cs > = «(R) (s'X ) 

-00 -oo 



C2.3. 1) 



Aa we know that UCR) is a radially aynuaetrie function and 
that the FT of a radially aymmetric function ia alao radially 
symmetric, the function K Cr cos r sin is radially 
symmetric. Therefore, without loss of generality, we may take 
^ = 0. Thus (2.3.3) reduces to. 
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n 

K C r > = q Cr cos © > d S 
0 

= CWr (r) 

i.e. , K is inverse FT of U. 

Writing K^(r) = (U^ we have, 

K^Cr) = K CAr) 

Hence we have, 

Z Cx 5 = C P U Cx) 

= C Z K (X) 

= C Z * K 5 Cx) 

and with K replaced by K^, 

(x) = CZ * ) Cx) ...C2.3.4) 

Thus, = { Z^Cx) : X e K* }, the ensemble of all 

reconstructed images using CBP algorithm is also a spatial 
stochastic process. 

2.4, VARIOUS SPACES 

5£ 

For a real valued function g on K , we define, 

g^ ^ Cx) = eCx^ + r cos x^ + r sin ^) , 

2 . 

where re [0, cb) , ^ e [0, ft] and x e R . 
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Let 


2n 




1 

2n 


g (x) d^, 

r , ip 


provided g ^CS) integrable with respect to tp . 
r f ^ 


DEFINITIOSF 2.4.1 The k-th radial derivative o£ function 


^ Cx) is given by 


“r,5 * ^ 


V r > 0 , X « K , and 


k = 0,1,2,.. 


Uhenever interchange of the order of integration and 
iio^ 

differer^ion is permissible. 


2n 


D -T g (x) = -TT- 
c ,ip ^ 2n 


JL Cx) dip. 


dir 


k ^r,iP 


0 


Let p be a positive integer. In order to further develop 
the theory we define the following spaces:— 


DEFINITION 2.4.2 The Space is the class of 

and continuous real valued functions g on K , 

g tCx) is 4p times continuously differentiable 

r,^ - 

Ic , 

r « [0, as), and D ^ gCx) are bounded and 


all bounded 

for which 
function of 
continuous 


47 


+ z 

functions of (r, x) in K x K for k=l , 2 , . . . , 4p . 

DEFINITION 2.4.3 The Space is the class of all second 

order wide sense stationary spatial stochastic processes such 

s 

that (x) is a bounded and continuous function on K , with 
the norm 

11 11 0 "" 11 Be "" X 1 1 

Note that j| 11 0 ~ Z(x) = Cp3. a quantity independent 

of X. 


DEFINITION 2.4.4 

order wide sense 
such that C^Cx) 
defined as 


The Space is the class of all second 


stationary spatial stochastic processes, 
belongs to the class The norm in is 



4p 



I D 


4p 

r, ^ 


c 


DEFINITION 2.4.5 The Space ^2p class of all radially 

symmetric functions K(r), satisfying: 


ca 


00 


CD 


r-^KCr) 1 dr = 


r^ 1 KCr) I dr < as 


V j = 1,2, . . . ,4p+l 
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OD 


(ii) 2n 


r KCr) dr = 1 


ca 


Ciii) 


KCr) dr = 0 


¥ j == 1, Cp-1) . 


Civ) 


m 

KCr) dr = 0. 


0 

and K Cx) = K C Ixj] ) « ^4^- 

DEFIWITIQM 2.4.6 The Peetre’s functional KpCt^^, ^ > is 

defined as follows: 


Kp C ^ > = inf ^ I 11 ^ -^0 

where 3^^ ^ ^n» ^ ^ * . t e K . 

0 0 4p 4p 


DEFINITION 2.4.7 For 0 < /5 < 4p, the Space 
4p 


is a space 


intermediate between spaces and tS. , which 

0 4p 

"interpolation space”. It is the class of all 
order wide sense stationary spatial stochastic 


we call an 
such second 
processes ^ 
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belonging to which can be decomposed in two second order 

wide sense stationary spatial stochastic processes, such that 
one of them (say, belongs to Sq and the other (say, ^ 4 p) 

belongs to and such that t ^ Kp C t^^, SC ') is bounded. 

rtr 

REMARK : Note that is not empty unless 1?^ is empty, 

because any spatial stochastic process 3C in Sq can be 
trivially decomposed as ^ ^ + 0 . where 0 is the zero-valued 

spatial stochastic process which obviously belongs to 

2.5 :^ME PRELIMINARY RESULTS 

LEMMA 2.5.1 If ^ is second order wide sense stationary 
spatial stochastic process^ then SC^ is also second order wide 
sense stationary spatial process. 

PROOF :— Using definition of as in (2.3.4) 

E C (X). (y) > 



Z (u) Z (V) C X 


u } C y ~ Y 5 du dv 




2 


= E 


Z (u) Z (u + w) K^ C X - u 5 X ~ H ~ ^ 


DR^E 


2 . 
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EZ (u) Z (u + w) K^Cx-u>K^CY~B~b?5 dudw 


IR^K 


2 


Cj| W i 5 c X - u ) Ka C y - U - W ) du dw 


JR^K 


2 


^11 y I 5 Ka C X - Y “Y 5 ^A c- V - w ) dv dw 


[e“ie 


2 


hCY5^j^C?“Y~Y5 dv 


= C X - y (say) 
where 


hC V ) 


C I V I > ^A C Y ■*■ If ) ‘^w- 


{R 


2 


For atat ionarity it ia neceaaary to ahow that 

^ ^ ^ ^ II 5 I > ■ 


First, we will show hCY5"^C|Yi5 
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Let V = C V cos B, v sin v=|jv|,S = tan 

and w = C w cos w sin w=||w|,^ = tan *C ) . 

Using these notations, h C v ) can be written as 

2n £0 

C w 

0 0 

K^Cwcos ^ + V cos &f wsin 4^ + v sin & >w dw d^ 

In m 

C w ) . 

K^Cv + wcos (< 5 ^ - ©), V + wsin - S) ]>w dw d4> 

= hCv)=hC|Y|>. 

Similarly, 

^AC?S"y^=j'5^C|lY(|)K^Cx-y-y) dv 

Hence E C Cx). (y) ) depends only on the distance 

between x and y, i.e., is second order wide sense 

stationary spatial stochastic process. 

CENT'^'t ' 

I ' 

' 11 ^ 


5sr 

•'•“■ H I M 



0 0 
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DEFINITION 2.5.1 The difference process ^ is defined as 

= { Z (x) - Cx) : X s k"" > , 
for Z e iS;, and Z^ e ,2:^. 

LEMHA 2.5.2 If -2^ is second order wide sense stationary- 
spatial stochastic process, then ^ is also second order 

wide sense stationary spatial stochastic process. 

PROOF :- E C Z (x) Z^ (x + u) ) 


ZCxJZ Cx + u-v5K^Cy^^Y 


K 


ECZ(x)ZCx+u-v>>K^Cv) dv 


K 


C I « - Y II > C V ) dv 


Let u = C r cos & , r sin B >, and v = C v cos v sin ^ > . 
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Then 

E C Z (x) C? + u) 3 

2n as 

= Cr- ~ V cos r - V sin - S) 3* 

Cv 3 V dv 

= gj^ ([ II u I 3» (say), a function of r only. 

Similarly, 

E C Z^ Cx) Z Cx + u) 3 = £2 € II B I! ^ ' (say). 

By Lemma 2.5.1, is second order wide sense stationary 

spatial stochastic process. Hence 

E C CZ - Z^ 3 Cx) CZ - Z^ 3 Cx + u) 3 

=E C Z Cx) Z (X + u) 3 + E C Z^ Cx) Z^ Cx + u) 3 

- E C Z Cx) Z^ Cx + u) 3 - E C Z^ Cx) Z Cx + u) 3 

= C II s u + C 1 a li 3 - Si C I « 1 > 

A 

- £2 C II B II > V X « , u e K® , 

which is a function of | u | only. 

Hence S - is second order wide sense stationary spatial 
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stochastic process. 

REMARK : Since is a second order wide sense stationary 

spatial stochastic process, the ESE 

Cx, =e[<Z-Z,K.^C2-Z^) Cx,] 

is independent of point x « ER . As it depends on R^ or A, it 

2 

will be denoted by E^ C 5 now onwards. 
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CHAPTER 3 


DIRECT THEOREMS 


3.1 INTRODUCTION 


In chapter 2 three types of error^i in reconstruction of 
an object image by CBP algorithm were defined. In this thesis 
our main emphasis will be on the theoretical estimates of the 
ESE, E® C ^ 

In this chapter we will discuss some direct estimates of 
the ESE for the object stochastic processes belonging to the 
classes defined in section 2.4. It is assumed that window 
functions used in the reconstruction by CBP^. belong to the 
class . 

2p 

The results established here use the derivatives of the 
covariance kernel of the object stochastic process and a 
modulus of continuity, which is appropriate for the analysis 
of ESE in CBP, i.e. , 


iii 


[k] 

C ; 6 > 


= sup 
0<|hl <S 

X « 




2 


(3.1.1) 


where h = h^ 


h , and h. , h « K 
2 1 2 
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3. a LEMMAS 


Some lemmas are given het-e which will be needed fot- 
proving the direct theorems. 


LEMMA 3.2.1 


[k] 


lim C ; <5 > 

S—i.0 ^ ^ 


0 if and only if 


lim - C„ (x) = T C„ (x) for every x « 0?^ . 

[ h|-.0 ^ 4 > ^ - 

[ k] 

This result is evident from the definition of <<> C C^. ; <5 > . 

Similarly another easy consequence of definition is 


LEMMA 3.2.2 6^ < 6^ implies that 

[k] [k] 

for k = 0,1,2,.... 

LEMMA 3.2.3 For any X and S > 0 , we have 


[k] 




Ek] 

XS y < Cl + X) 6J^ ( 


5 ) 


PROOF 

Without loss of generality it suffices to prove the Lemma 
for the case h > h^, h , h « CE^ . Ue have from lemma 3.2.2 

1 a ji St 



ts^ 
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[k] 


[kl 


6ic C ; X6 5 < c ; (1 + [X3) ^ > 


Now, W6 have, with h = h^ ~ ^2 


[k] 

C ; Cl + [k]) 6 ) 


sup 

0 <h< Cl+CX]) & 

X <s D?* 


“h . ^ 's> 

i 


°h . ^ 
2 


sup 

0 <h< C1+[X]) & 

, t ->2 

X « ER 


K 2 Cx) 

h +h, 4> ~ 

2 


2 


max 4 sup 

1 <}< l + CX]'^Cj-l) <5<h<j 6 

m2 

X « IK 


^h .h, 

2 


2 


»i) 


max -1 sup 

1 <j< 1 + [X]^ 0<h< S 

X s 




C w • . N,- rC_ Cx) 
h^ + Cj-l)h,^ :Z 


1 } 


l + [X3 




sup 

0<h< S 

X s E* 


Vih 'i’ 


D, . , . ^ Cx) ^ 


'h^+Cj-l)h, ^ 


[k] 

C 1 + [k] > “c ^ ; *5 ) 


[k] 

< C 1 + k > C ; <5 > 
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LEMMA 3.2.4 The integral-type remainder in Taylor’s 

expansion is given by 





r. 


k 

fs’ - 2 iT “o, f 

j = 0 


1 

Ck-1) ! 





Cx) 



- c 


} 


(x) Cr 


s ) ds . 


3.3 THE 6i - ESTIMATE 
e 


For a second order wide sense spatial stationary process 

s: , the following theorem relates the ESE E* C > to the 

smoothness of the covariance kernel C^. This thoerem gives an 

estimate of ESE for the cross-sections belonging to those 

processes which have derivatives of covariance kernel with 

[k] 

property i.e. C ; <5 ) > 0 as «5 >■ 0. 


THEOREM: 


Let the CBP window K e K- . Then for k = 0 , 1 4p-l ) 

2p 


. ®k t*'' 1 

4 ^ 7E- “= « ■ A > 

A 

where B is a positive constant independent of and A; 

k 
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provided 


PROOF: 


Ek] 

;< 5 > *.0 (<5 ^ 0 ) 

= E C - Z >^Cx) 

2fi 2n m m 

C ZCx^- I cos S, x^- I sin - Z(x^,x^) >. 

0 0 0 0 

C ZCx^- cos x^- I sin - ZCx^ .x^^ ) > . 
r K C r ) s K C s ) dr ds dS d^ 

Zn Zn c£i ca 

E C Z ( x^- cos e, x^- sin & ). 

0000 -7^ ® ^ s... 

Z C X - -^ cos X - -j sin <jiS> ) 

i A Si A 

- Z ( X - cos S, X - V sin S ). Z (x , x ) 
^lA ’aA 

- Z (x - cos X - f cos ^ ).Z (X , X ) 

1 A 2 A A 2 

+ Z (X^. X^) Z Cx^, 2C^) >. 
r K(r) s K(s) dr ds dS d^ 

Zn Zn ca CO 

Cp, C (r-* + s* “ 2rs cos / A 5 

0 0 0 0 




m ^ 
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- C^. C r/A > - C S/A > + C 0 } j 

r K(r) B 1 ( 3 ) dr ds dB d# 


2fi CO 


= 2n 


^ 2 S 

C Cr + s'" 


0 0 0 


2rs cos (9)*'"* / A > 


- C i^/A > - C S/A ) + C 0 > j 
r KCr) s K(s) dr ds d^ 


Using Taylor’s expansion given in Lemiaa 3.2.4, for 


C -/A > = .,C 0 >, 


3^, -/A^ 
2n (£! i£> 


c ;« ) = 2 n 


j^JCl/A'^ji 5 (r® + s^ - 2rs cos & 


0 0 0 




D 


j _ 



r K(r) s K(s) dr da dS 


2n CO <u p/k 
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r/A 


‘^2 f 2 ) ■'*= 


“S, i =2 '2) >.ci - t)‘'-ldt 


0 

s/A 


’■K.i C2) 


■>1^ t2) 5 


k-1 1 

t) dt 

J 


r K (r) s K(s) dr ds dS , 


2 S S£ 

where p = r + s - 2rs cos &. 


Since Cx ) = C f x | 5 is an even function of |[ x j| , 


D- - (x) will be zero for k odd, and since K « tK- 

0 , 4? ’ 2p 


e!c^) = 1 eJc^) 1 


2n CD 00 p/A 


1 2n 



r 

1 Ck-l)l 



1 


0 0 0 


C®! 4 . '23 - <, J (5) ). 


C f - t dt 


r/A 


C2) - i>S, ^ (2) > 


0 
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S/A 


'9) 


CO) ) c-| - 


2n tib m 




Zn 


Ck-D! 


I “1 


r K (r) s KCs) dr ds d6 


^ '9) - Car <-<>) |. 


0 0 0 0 


! f “ t ^ dt 


r/A 


I “1 i Ca (0) - D^, ^ CO) I I ^ - t 


S/A 


I “t.^ '9) - 


<2) I l| - t|'‘-^dt 


r K (r) KCs) dr ds 60 ...(3.3.1) 


Now, the integrand in (3.3.1) is 


C‘'3 p 

^ “c ( C^, ; e > 


P/A 

r 


k+l r 

t dt + C ^ ) 


0 
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t ^ dt + CC^ = 5 ^ I 


1 f o ^fi r r k 

|[ -c ^ I > ^ f > " -c ^ ^ A > ^ A > 


^ ^ k 1 

c ^ ' A ^ A j 


^ [ C 1 + P 5 1 ^ r ) . < 1 * s > =‘‘ ; . 


"= < '^;* = S > • 


(by Lemma 3.2.3) 


Hence . 


9 I "■ J ■! 

E, c ^ > < — tj c ; 4- >. 

A ,,»k c -=£ A 

k i A 


2jt CD CD 


C (1 + p); 2 .^ + (1 + r) r^ + (1 + s) s^>. 


0 0 0 


r j K(r) I s j K(s) | drdsdsj 


^ Zn 1 > 

" k “c ^ ’ A ^ ‘ 


kJ A 


2 n 05 OD 


-, 2.2 - 

C (r +s - 2rs cos 6) 


0 0 0 
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Ck+l)/2 

+Cr +s - 2rs cos &) 


r I K(r) j s } KCs) j dr ds da 


00 CO 


+ 2fi 


k+ 1 k 

C r + r > r j KCr) | s j KCs) j dr ds 


0 0 


00 00 


+ 2n 


Cs^^^+ s^ r j K(r) j s j KCs) j dr ds J 


0 0 


where 



[k] 

C ; 



» 


€£} m 

Ir 4- i k 

C S + s > r j KCr) | s | KCs) j dr ds 

0 0 


= 2 


C4n) ' 
k ! 


2n 00 00 


. C4«) 
kl 


C Cr*+s^- 2rs cos a)^'^^ 


0 0 0 


^ 2 Cki-l)/2 

+Cr +s - 2rs cos a) >. 


r j KCr) I s I KCs) | dr ds da 
Note that since K <e • 

Thus the existence of higher order derivatives of covariance 
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kernel with io property implies less ESE, because o£ A 
c 

occuring in the denominator. 

3.4 ASYMPTOTIC FORMULA 

The following theorem in this section gives an estimate 
of ESE for those spatial stationary object processes SL 
possessing a covariance kernel which is very smooth in the 
neighbourhood of origin. 

THEOREM 3.4.1 


4p 

Let K « K- , ^ and suppose that D r C~ (x) exists 

2p 0 r , ^ - 

at X = 0, for r « [0, as ). Then following asymptotic formula 
for ESE holds: 


E^C 


r jlny 1 
I C4p)! J 




•4p 


+ o C 




as (A ► as) , 


PROOF : - 

Since t (0) exists for r s [0, as ). The Taylor 

expansion of covariance kernel ^ ^ (0) about r = 0 can be 
written as follows: 

4p 

fr “o. f 

(3.4.1) 




CO) = CO) + J 

j =1 
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where 


Ap 

1 

j=C 


■"a = « ^a.r, ^ 'S) - J ^ dJ_ ^ C^ CO) ) / r 


4p 

(3.4.2). 


Ue note that. 


1 iai . 
r-+0 ^ 3 ^ 


Cr) 


0 


(3.4.3) 


and 


liffl 

r-+aj “ar 


Cr) 


exists . 


(3.4.4). 


Since C - (.) is a continuous function of r, h~(r) is also 
continuous. Hence by continuity of h^ (r), (3.4.3) and 
(3.4.4), the function h^ (r) is bounded. Let H be the bound 
of (r) , i . e. , 

I h^ (r) I < M, for r « [0, co ). ...(3.4.5) 


Now, 

E® ) = E C - Z (X) 

2fi 2n m m 

C Z(x^- I cos &, x^- I sin S) - Z(x^,x^) 

0 0 0 0 

C Z(x^- cos x^- I sin <^) - Z(x^ ,x^) 5 
r K ( r ) s K ( s ) dr ds dS d^ 


= E 
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where 


4p 


Cr) = C ^ ‘5) - 2 -7T ®0. CS) > / r 

j=£) 


4p^ 

. . .(3.4.2). 


Ue note that, 


lim 
r— i-0 


(r) 


0 


(3.4.3) 


and 


1 iui . 
r — *&} X 


(r) 


exists . 


(3.4.4). 


Since ^ ^ (.) is a continuous function of r, l^(r) is also 
continuous. Hence by continuity of (r), (3.4.3) and 

(3.4.4), the function h^ (r) is bounded. Let M be the bound 
of h^. (r) , i . e. , 

I h^ (r)|<M,forr4E[0,£o). ...(3.4.5) 

Now, 

e“ > = E C - Z (X) 

2nr 271 CD GO 

C Z(x - I cos S, X - f sin <9 ) - Z(x ,x ) )• 

0 0 0 0 

C Z(x - ^ cos X - f sin - Z(x ,x ) 5 
'• iA zA i ’ 2 ^ 



r K(r) s K(s) dr ds d8 d4> 



Since 

gives 
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Zn 2n (it CO 

E C Z ( x^- cos &, x^- sin ^ )_. 

iJ mJ ■. 

0 0 0 0 

z C x^- ^ cos x^- •! sin ^ ) 

- Z ( x^- ^ cos e, x^--^ sin e ) Z (x^, x^ ) 

- Z (x^- ^ cos x^- j cos ) Z Cx^ , x^ ) 

+ z (x^, x^) z (x^, x^) 

r K(r) a K(s) dr ds d6 d(^ 

m <30 Zn Zn 

C + s* - 2rs cos C& - / A > 

0 0 0 0 

- C r/A > - C S/A > + C 0 ) j ^ 
r KCr) s KCs) dr ds d^ d<?!> 

SS <£ and K e application of Fubini’s theorm 

<30 CO Zn Zn 

CCr* + s* - 2rs cos (& - 4^'if''^/ A 5 

0 0 0 0 

- C r/A ) - C S/A ) + C 0 ) j ^ 
r K(r) s KCs) dr ds dS d<^ 
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Zn 


In 


[ C + s* - 2rs cos / k ") 


0 0 0 


- C r/A > - C S/A ) + C 0 ) 
r K (r) s K (s) dr ds d& 


oi 60 2n 


2n 


^ ":*,t/A. S (2> 


0 0 0 


,A CO) 

^,r/A, 


- ‘'*,s/A, i, 'e> * t2> 

r K(r) s K(s) dr ds , 


. _^2 2 2 _ . 
where t = r + s ~2rs cos & 


Therefore, Taylor's expansion of - (0) about r = 0 gives 

, r 


<3b <£t 2n 




2n 


4p 


[2 


j!A 


j C s*-2rs cos ~ ), 


j = l 


0 0 0 


“j. ? (2> * •'a « I 5 

- Cr/A)^^ ^ I > -Cs/A)^^ ^ S ^ ]. 

r KCr) s KCs) dB dr ds 

Since as A ► 60 , r/A, s/A, t/A — > 0, using 3.4.3 and (3.4.5) 


we get that 
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r 4p+l 4n+l 

(i) [ (^ ) (r) K(r) j is bounded by M (r/A) | K(r)j 

an integrable function. 

(ii) Since C^) >■ 0 as A > oo, 

Cr/A)^^"^^ C^) K(r) — + 0 as A 

(i) and (ii) imply the applicability of Lebesgue dominated 
convergence theorem. Hence, 

00 

Cr/A)^*^"^^ h^ C|) KCr) dr i- 0 as A ^ . 

0 


Thus , 

00 OO 00 

[ Ct/A)^^ ^ > ). 

0 0 0 


r KCr) s K(s) 60 dr ds ^ 0 


as A > m 


Hence , 


00 00 2n 


E‘ ) = 2»t 


4p 


(n 

j=i 


2 2 


lA 


j CCr*’+ s'-2rs cos r"^ - s'*>. 


0 0 0 


:o3] 


K 2 CO) I r KCr) s KCs) 60 dr ds 
0 , <P ^ - 


+ oC ) 


(A + 00 ) 
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Recall that C^. ^ ^ (x) is an even function of r, therefore 
j CO) = 0 for j odd, and since K « K„ we have 

0 , # ^ - 2p 


(3Ci (£> Zn 


2n 


4p 


( 27^1 

j = l 


2 2 


j C Cr^H- s‘‘-2rs cos r^ - ). 


0 0 0 


Dn 2 (0) 1 r KCr) s K(s) m dr ds 

u , ^ -® - J 


CO oa Zn 


= Zn 


2p 


0 0 0 


(2 

j = l 


A® '* ( 2 j ) ! 


C r®’’+ s“^+ C-2rs cos 


i=i 


2 j —2 i -2 2.i.^ v-, 

s ^ + Cr + s ) C-2rs cos d)'’ j 


Dn"' 2 Co- CO) j r K(r) s Kfs) dr ds 


)2j 

^0, 


00 CO 2fT 


= Zn 


( “o! i '^2 ‘2) / a’^C^p)! ] ( 


> Crs)^^ + 


0 0 0 


p-1 

C2rs cos + J C2cos «) 

i=l 


(rs)^*^ 


]. 


r KCr) s KCs) dB dr ds 

CO CO 


2" ( “S! i <5)/A^'’C<P)i ] 


Crs)^^^^ KCr) KCs) dr ds . 


0 0 
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2n 
' P 

2 <2i ^ > ^2cos S)2p-2i da 

- i = 0 


[ ^ ] ( can)""/ C4p)i 


[2^a? 


= 0 


> c^j ) c ) ] [ 


]• 


m 


kCr) dr 


[ ”l ^2p / A ^ ] ^a. ^9)' 


where 


= [ 


= I can)'/ C4p) 


■] 


C > I and A- 
P J ap 


<sa 


K Cr) dr 


Hence, 


E*C^)=0CA^^) 


as A »■ CO 


3.S THE I C^ ESTIMATE 


THEOREM 3.5.1 


I£ ^ e te. and K e K , then 
4p 2p 


E® C ^ ^ C H2 / > I C 


a «4p ’ 
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whece 

as CO 2fT 

Cr^+ s^- 2rs cos B rs jK(r)j (KCs)} 60 dr ds 

0 0 0 

mod 

+ 4n 

0 0 

PROOF 

For any x « K*, 

> = E C - Z (X) 

2n 2n w i3b 

C Z(x^- I cos &, x^- I sin & ') -ZCx^.x^j^) >. 

0 0 0 0 

C ZCx^- ^ cos x^- ^ sin -ZCx^.x^) >. 
r KCr) s K(s) dr ds dS d^S- 

2n 2fi m <3b 

E c z c x^- cos e, x^- sin 9 ), 

0 0 0 0 

Z ( x^- cos x^- ^ sin ^ ) 

- Z ( x^- J cos e, x^- ^ sin 6 ) Z (x^ , x^ ) 



r4p+l ^ I (r)| jK (s)| dr ds 

.j 




m ^ 


7 3 


■ 2 ( x^- cos 4'. ^ sin ^ ) Z ) 

+ Z (x^. x^) Z (x^. x^) ). 
r KCr) a K(s) dr ds d^ d^ 

2n 2n m (£t 

- r r r 

[ g 4 

C Cr + s - 2rs cos (S* - ^)) /A > 

0 0 0 0 

- C^- C r/A 5 - C S/A > + C 0 ) j . 

r KCr) s K(s) dr ds dS d^ 


m €io 2n 


= 2n 


^ . 2 2 

C (r + S 


2rs cos / A ) 


0 0 0 


- C r/A > - C S/A ) + C 0 > j . 
r K (r) s K (s) dr ds da 


Since ^ e implies that s Taylor’s expansion with 

t^ = r^ + - 2rs cos a gives, for some s (0, 1), i = l,2,3 

CO oo 2n 


> = 2nr 


P 4p- 1 


(17^ 

j = l 


. 2 St ^ stX2 j ,.j. 

j CCr + s -2rs cos a) - r - s O 


0 0 0 


^ ^ l/A^^C4p)i > C t^^. 


2 ( 0 ) 
ol^t/A,^^^- ^ - 


« 


r/A,^- 
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r K(r) s K(s) d& dr ds 
ca Qd 2n 


2n 


= c 


A^^C4p) ! 


J J J 

0 0 0 


r ^ 2 2 _ >-, ^2p 

(r + s - 2rs cos ©) , 


:* C_C0) - r (0) 

ot^t/A,^ oi^r/A,# j£ - ' 

- .J. T C_ CO) 1 

a s/A,^ ^ J . 

r K(r) s K(s) da dr da 


because D;i t (0) = 0 for j odd and K« !K^ . Hence 


Ci^>= EtC^> 
A 


2n 


€£t (£f 2fl 


‘a^^C4p ) I 


[ Cr*+ s^- 2rs cos j rs . 


J J 

0 0 0 


l°»'’t/A,^ '^3^2^! 


00 m 


+ 2nr 


0 0 


00 00 


+ 2n 


4p+l j 

& : 
m 

1 KCr) 

1 ^ 

,4p+l j 

jj4p 

at 

a 


C~ (0) i. 


2 '=3 tS) I, 


0 0 
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KCs) I r I K(r) j ds dr J 


00 as 2 n 


2n 


^ c 


A^^C4p)! 


•> i i 0 


Cr^+s^-2rscos 


0 0 0 


rs [ KCr) I 1 KCs) j dS dr ds 


W 00 


^ 4fi 


4p+l 


s j KCr) [ I KCs) j dr da j 


0 0 


C H2 / A**’ 5 II J lo 


Since K « ^2p’ 


we have M 2 and thus 


e' < a > s c tij / a’p ) 1 ^ io 

^ < «2 1 «4p 

Note that, this theorem gives an estimate o£ ESE for those 
stochastic processes which belong to the class • 
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CHAPTER 4 


INVERSE THEOREM 


4.1 INTRODUCTION 

In chapter 3 we have seen some estimates o£ the ESE 

2 

M > for CBP , whenever the object process M is known to 
have certain theoretical properties. For exampley either the 
process belongs to some class of stochastic processes as 
defined in chapter 2, or it has soiae other properties such as 
smoothness of its covariance function and the conditions 
related to the modulus of continuity. 

In this chapter we solve the inverse problem^ i.e.y if 

the order of ESE is known then to which class does the object 

process M belong? The result in this chapter states that the 

ESE is of order 0 ( A if and only if ^ 

interpolation class, provided the CBT window K e ^2p' This 

theorem guarntees the error order to be A ^ if the object 

process can be represented as a member of the int erpolat ion 

class ^ . Note that as ft increases the covariance function 

4p 

of the object process becomes smoother- 
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4.2 LEMMAS 

Ti'i6 £ol lowing Leaiflias air-e f-aqiiif-ad in the pt-oof' of the 
main theorem. 


LEHMA 4.2.1 

Let -S: s K « fj [ 0 , co ) and 

4t> 1 

^ KCx) « L [ 0, ® ) and be bounded. Then 

PROOF:- 



wide 


For proving that 2^, « 

A 4p 

€ ^4p> because Lemma 2.5.1 
sense stationary process. 


we need only to 

states that is 
A 


show that 
second order 



E C Cx). Cx + y) ) 


Z (u) Z Cv) (X 


u) <x + y 


v) du dv 




E C Z Cw + V) Z (v) > Cx - w - V). 


(x + Y - v) dw dv 





7a 


^ 5 K, Cx - w - v) K. (x + y - v) dw dv 


2 2 

K K 


j C5f) C-w - u) (y~ y) dwdu 
IR® DR® 


J Cu) (Y - u) du 
K® 


where 







C^ Cw) (w + u) dw 


m 


2 . 


NOWy 


I CU)| = I 


(w) (w + u) dw 


DR 


= I 


(w /A - u) K (w) dw j 


IR 


2 




.1) 
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< 


I Cw /A " u) j |K (w) I dw 


ER 


2 


I I, 


I K (w) j dw 


0? 


Hence (u) is bounded in K* . 

2 

Let u, V « ER , 


1 ig (say). 


j (u) - Cy) I 


= I 




w /A ) - Cv - w /A ) ) K (w) dw j 


ER 


^ I Cu - w /A ) - Cv - w / A ) j 1 K Cw) j dw 

IR 

Since integrand is bounded by 2 | [[p j K (w) an 

integrable function. Applying bounded convergence theorem, 


Hv I “a «S) - h* Cv) I 

- (y - w /A )| . 

dw = 0 , 


lim 
u — 


I C Cu - w /A ) 


ER 


K (w)j 



so 


Since ^ C^. is continuous in . This implies that is 

a continuous and bounded function in 
From (4.2.1) 




(u) K, (y - u) du 


£R 


Cy - U /A ) K (u) du 


[R 


Hence, on similar lines, C^. is bounded by H jj „ and 


continuous in 0? 


Since is bounded and continuous in K and K (x) e fl, , and 

4 1) 1 tf 

JJ F K (x) -eL [0,oo), D t C^. (y) exists and is bounded 
r ,<p ~ r,<p - 

“I* ^ 

and continuous function of (r, y ) e K x K , and is given by 


D 


r, i> 


(u) A 


2 S 


Sr 


k < "^Ar. ^ ( A Cy - u) ) > du 


Since K (x) = K ( I x | ) « [0, ) and 

- K (x) « [0, 00 ) the integrability of ^ C?) 

follows from interpolation inequalities [see Goldberg(1966) ] . 
Hence, e ^4p that e ^4p' 
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LEMHA 4.2.2 


Let ^ K (x) € L [ 0 , 


ac} 


then there exists 


constant such that 


e| C a > < (2 113)^ 



I 


0 


PROOF:- 


> = E C - Z )^Cx) 


= E 


C Z Cx - u /A ) - Z (x) > c z (X - V /A ) 


2 2 

E E 


- Z (x) y K (u) K (v) du dv 


= J j C C Cu - V)/ A > - (u /A ) - Cy /A ) 
E® E® 


+ CO) ) K (u) K (y) du dy 


Hence , 


E^c ^ 


I e!c ^ ) I 


[ ! C (u - y) /A > I + I Cs /A) 1 + 


■^2 „2 
E E 


I C, 


I C 9 > I ^ Cu) K Cy) du dy 


Cy /A) 


"¥ 


S2 


^ i io ( i ^ I ] 




C 2 M 3 > I II Q 


where 


M 3 = 


j K (u) j du 


ER 


LEMMA 4.2.3 

Let SL and ^ be such that the process ^ e . 


Then 


I 


Sf + ^ HQ 


in ^ 2 C i C, 




LEMMA 4-2.4 


Let « Sq, 


K € L^ [0 


00 


then there exists a constant 


M 3 < 00 such that 


wh^r^ «3 



in Leinina 4.2.2. 


LEMMA 4.2.5 


Let ^ K (x) « 


L^[0, 


00 ) , and D 


4p_ 

r,4 


KCx) -e L [0, ®) 


Then there exists a constant M^ < ao independent of A and ^ 



S3 


such that 

i ~ ^4 II io' 

PROOF :- 

(y) = E C (X) Cx + Y) 5 
A 


= E 


Z (u) Z (v) (? “ C? ■'■ ¥ “ Y) 


S 2 

[R IR 


:j J E c z (w + u) Z (v) ) Cx - w - V)- 

K® (x + y - V) dw dy 


Cx - w - y) (x + Y - y) dw dy 


ER 0 ? 


Cv) K, C- w 


uD Cy 


u) dw du 


2 2 

[R 0 ? 



IR 


(w) (w + u) Cy - u) dw du . 


c K,„ 7 c A (y - u) > > dw du 


SCKr') 


4p ' “Ar, # 


= A 


4p 


4p 


Cw) K^Cw + u) ^ K^CY - u) dw du 


ER® 05 ® 


Heiice . 




CR® 


I Cw) I I K, Cw + u) J. 


i ^ ^A I 


i C^ Iq Ms L, (say) 


this implies that 




l d‘p 

« r. 

- c 

i 

where 

L < CO 

and 

< <Xj. 

Uaing 

Lemma 

4.2.4 


1 

*0 

ti 

1 1 

k 

0 ^ 8 


4p 


!IO “3 


^ Ss:. 8 o 


- ^3 S 1 0 ^ ”3 ^ 8 Sk 8 0 

< H A^*^ I 80 ^ ”3 ^ ^ ^ because A > 1 
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Using Leaiffla 4.2.4 


IL.. = I lU + I r C 


»4p 


“0 

A 


r, 4 ^ 8 0 

4p 


^ «3 I lo " "3 I K.i C* «0 


= n c 

3 


4p 


4.3 MAIM THEOREM 


THEOREM: 


Let JS" e Sq and suppose that K « ^2p‘ 


C ^ ) = 0 CA~^) 


as A 




if and only if for 0 < /5 < 4p. 

PROOF :- 

First, let C ) =0 as A »■ ® . Then we 

ft 

to show that ^ « ^4p’ ^ decomposed as 

= ^0 ■" ^0 ^ ®0' ^4p ^ ^4p 

and K„ C , ^ > = 0 (t^) as t 0 

1: 

Let us define another ifeetre’s functional as: 


* 4t> 

Kp (t^P. ;ar ) 


"a^ ^ - X h * i i4p ^ 


Thus 


have 
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Kp ( a: ) < K* c :z y 


Cx) = C ^ >^Cx) (x) 


and hence, 


Cy) = C 


C;a; - 


E (Z - Z ) (x) (Z„) (x + S’) 


+ E CZ.),Cx) (Z - Z„), (X + y) 


^css: - sc ) ■*■ T*--^ 

k,r, 4 > B A,r,^ B'^A,r,<^ 




'£’ '^B^A.r, C! * !f> '‘■i* 


E C 2 g)^ (X) CZ - 2 g)^_^_ ^ CS + i) d.,! 


2 2 

Let Cx) be a measure on K such that ^ (K ) = 1, and 
Y (x) = Z C-x). Then, 


^ E CZ - Zg)^ (X) (Zb^A.p. ,6 (X * 
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2n 

E (Z - Cx) CZg)j^ j, ^ (5 + y) C?) d4> 

0 £R® 


In 


= E 


- 2b)a -i;? 


^^B^A,r, + I) C ?5 


K 



(Z - Zg)^ (X) ^ C? + Y) d^ (X) 


[R 


2 


= E (Z 


2n)A * CZt,), 2 C-y), 
B'^A (U ^ B'^A,r,«ii> ^ -■^' 


where * denotes convolution with respect to measure p. 


Also , 


2n 


1 

2n 


CZg)^ Cx) CZ 


^B’A.r, ^ 'i * y’ 


2it 

E (Zg)^ Cx) CZ - Zg)^^^^ ^ Cx + jr) d^i Cx) d^ 

0 IR* 
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2n 


= E 


C?) 2n 


+ y) C?) 


D? 


(Zg)^ Cx) (Z - Zg^A.r, ^ y) ^ (3?) 


ER 


E cz„), * C Y - Y„), - C-y) 

^ B-^A iJ ^ B'^A.r, ^ 


E (Y - Y„)j * cz„)j - c-y) 
B A ;u B A,r, ^ 


Hence using the fact that K -e 


2fi 


I D 


4p 

r. 


_1 

2n 


E (Z - Zg)^ Cx) (Zg^A^r, 1 


& 


4p 


ar 


^ E (Z - Zg)^.; (Yg)^^^^ ^ C-y) 




4 p 


0r 


^ (Z - Zg)^ * C^B^A.r, 
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a 


4p 


E (Z - 2„). * — CY„), - C-y) 


~ # ^^B^A (? + y) (?) 1 


0? 


E (Z - 


2g)j^ (?) ^ 5 ^ 


K 


C E CZ - Zg)^ (x) 


CR 


C E (Zg)^ Cx + y))^ dAi Cx3 




K 


C j E (2^)^ C? + y))^ 1 (x) 

(4.3. 


and 




1 ) 
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.4P 


‘'A.r, ^ 


^ T (x + y - u) du 


1 

u 

J • 


IR 


# 


4p 




IF ^B.r. ^ * X - X) “a <‘5! 


K 


^4p 

Z_ C'i) « — Z_ ■T 

B ^ 4p B,r, ^ 




(x + y - v) (v). 


2 2 
R [R 




K 




4p A,r, 4 > 


(x + ^ - u) du dv 


»4p 


Sr 


^ EC ZgCu) ^ C? + y - V) (V). 


g , jg* 

ER K 


? a4p 

* K 




4p Ar, ^ 


C A (x + y - u) > du dv 


S 


4p 


4p T 

B,r, <p 


Cx + y - u - v) K- (v) . 


2 2 
CR E 


= A 


4p 


^4p 2 ^ 1 Cx + y - u) du dv 

r , # A , r , # “ - -* - - 


^ C? + y " S “ Y) Cy) • 


B 


E E 


D 


- K, - (x + y - u) du dv 

r, ^ A,r, ip - 
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Hence , 


® C “r! i 'S * 


< A 


4p 


I ; =a, cx * y - s - y) I I It, 


B 


[R® K* 


- K - 

r, # "^A.r, 4> 


(x+y-u) I du 


< A 


II - C li 

11 r, ^ ito 

15 


I CD j . 


K® IR^ 


I ■>r!’ } h.T. ^ ts * y 


- u) 


£ A*" I C„ I 


'*B tip 


I K Cv) I dv 


K 


I 2 ^ 2 . ^ T (x + y - u) j 
• r, # A,r, - - - • 


R' 


- «3 ■- 't*'’ I <^2 l,p 

Jd 


CD I . 


dv 


j du dv 


du dv 


. C4.3.2) 


Now using Lemma 4.2.2 and C4.3.2) in (4.3.1), we have 
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In 


I ^ 

I V, 4> In 


E cz - (x) ^ (5 + y) j 




1/2 


ER 


C a"**" M 3 L I C^ diJ (x) 

B 


C A^» m2 I 

B 


IK' 


< ”3 '■ I ‘ 2 agl 4 p 


C A ’’ Mj 1 C^2 _ a^5 Ijp >^<'^■-11 '^atgll4p ^ ^ 


< 


^ m 5 11 C,„ „ . i„ + H. L 11 C. 


( 11 


3 » '"(Z - sjg) HO ■ “3 “ " '':afgH4p 


(4.3.3) 


Similarly, 


2n 


j - -^ 
« r, # 2fx 


E (Zb)a (S3 (2 - Zb’a.p, ,s (S + X) I 


0 
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C2 M 3 ) A^^ll C 


H3 . I cm^ 


B 


B 


”5 « 1 c,. 


io a JUp ^ 

B B 


where M, = 2 M. + Mr + L M^ . 

3 H J 3 


Let ^ = i ^ ^ . 

** rk 


Since as - SC e ^ and ;g. s 

u AO A 4p 


K* c a > i II _ a to + II 1^^ 


and since, 


^SC - ;a- = E cz - Z^) (X) (Z - Z^) (X + y) 

A 


< CE CZ - Z^)^ (X) CE CZ - Z^)^ Cx + y) 


= E CZ - Z^)" Cx) 


Ea C ^ . 


and E Z Cx) is independent of x 


i ^SC - ^^1^0 ^ ^ i ‘“Sf^tUp 


< -E^ C ^ C C I C^ _ ^^i,Q . H 


.4p 


4 P 


^ I '“;af I ^ 

B B 


Finally using Berens-Lorentz Lemma, we have 


K* ct’”. a 5 = C I . a llo * I '^a.hp > 
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= c 0 ca“^) + ctA)^*" i ^ > 

= 0 Ct^) as t . 0. 

Hence M ^ . 

4p 


Conversely, let ^ This implies that 


^ = ^0 ■*■ ^4p' ^0 ^ ^0 ’ ^4p ® ^4p* 


C_ I < H, t 
•^n D 6 


ft 


and 


... (4.3.5) 


<'a,Jl4p ^ "7 


ft- ip 


Hence , 


E' C ^ > 


E CZ^ - Z)" (X) 


E (Z^ - Z^ + Z, - Z. ) (x) 

°A ° ^^A 


^ 2 E (Zq^- Zq ) (x) 4 E (Z^p^- Z^^ ) (X) 


Now, using Theorem of section 3.5 and Lemma 4.2.2, it follows 
that 

c * ) s 2 C C2 H3)" I €3,^10 + C H^/ A*”) 1 > 

Uith t = (l/A) in (4.3.5), we get 
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CHAPTER 5 


SIMULATION RESULTS 


5.1 INTRODUCTION 


Ift chapters 3 and 4, we have discussed some theoretical 
estimates of ESE in continuous convolution backpro jeetion 
method. But in practice the discrete version of CBP algorithm 
is applied on digitized images. In this chapter we consider 
this practical case of discrete implementation of CBP 
algorithm and give some estimates of ESE which incorporate 
the errors due to discretization and interpolation in 
addition to the inherent error. Thus, we shall estimate the 
error in reconstruction given by 

i j 

where, MP is the total number of pixels used in the summation 
and Z (respectively Z (x.j) ) denotes the value of 
function ( respectively reconstructed function ) at i,j-th 
pixel . 

Ue recall from Sec. 2.5 that for a second order wide 
sense stationary spatial stochastic process X, with 
covariance kernel C^, the ESE is given by 
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e! C 5 = E CZ - Z 


( Ml H - Y |) - Cl U I ) - c I vjl ) 


z z 

K D? 


■*■ C0)^KC 

(5.1.1) 


Let I be any set o£ indices and suppose that the values of 
covariance kernel are known at points d^ , i «e I . Then from 
(5.1.1) it follows that can be approximated by 


i€l 


. . . (5.1.2) 


where the co-efficients depend upon the CBP window used. 


Since by (2.3.3) 


we have. 


Z (X) = Z * K (X), 


(u) = E C Z (X) Z (x + u) ) 


Z (v) Z (w) K (x - v) K (x + u - w) dv dw 




= E 


Z(w + y) Z(w) K(x - w - y) KCx + u - w) dy dw 


K^CR® 
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where K_Cu) = K (-u). This implies that C~ is a linear 

48 & 

function of and hence we have a further approximation of 
ESE given by 

^ ^ = 2 ^i ...C5.1.3) 

iel 

where the co-efficients again depend upon the CBP window. 

If the values of the covariance kernel Cd.) are known 
for i « I, (5-1-2) suggests regressing these on the values of 
E^ C af > for a sample of object functions and obtaining the 
estimates of ot^’s. Using these estimated values of ct^ and the 
relation (5.1.2) we can predict the reconstruction error for 
the future object functions for which the values C (d.) are 
known for i « I. Even, if these values of covariance kernel 
are not known, they can be estimated and the estimated values 
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can be used to determine /5.’s in (5.1.3). Finally (5.1.3) can 
then be used to predict the reconstruction error for future 
object functions. 

The impl emetation of the above techniques of predicting 
the reconstruction errors in CBP algorithm is illustrated by 
carrying out a detailed simulation study. This simulation 
study also shows the validity of the u> - estimate of ESE. 


5.2 SIMULATION OF RANDOM IMAGES 


For generating a typical 64x64 stationary image the 
following method was employed : 

a) A 72x72 array of independent pseudo random numbers 
were generated from a uniform distribution on the interval 
[0,31]- These numbers when rounded to the nearest Integers 
are called the codes (CD) and represent the pixel intensities 
in 32 gray levels. 

b) For each of the inner 64x64 pixels, a circle with the 
radius of four unit pixel length centred on the pixel is 
taken and the weighted average of pixel intensities inside 
the circle is obtained and substituted as the value of object 
function in that pixel. This gives rise to a 64x64 stationary 
image with a covariance kernel depending on the weight 
function. 

In our study, we used the following eight weight 



105 


funetiona : 

(i) exp C-r), (ii) exp (iii) exp 

Civ) exp C-r^), (v) ( 1 + 10 r (vi) (1+10 

(vii) (1+10 r^) and (viii) ( 1+ 10 
Height function (i) was used for generating 70 

independent sample object images, and eight sample images 
were generated by making use of each of the other seven 

weight functions. Henceforth we shall refer to these sets of 
images as samples. Thus sample 1 consists of 70 images 
corresponding to weight function (i), and samples 2-8 
corresponding to weight functions (ii)- (viii) respectively 
consist of eight images each. The projection data for 
these images were then calculated on computer as follows : 

For each line ( 1 = 1 L, and k = 1, K, the 

projection data is given as integral 


T 


pCSj^, Sj^) 


z (x cos + t sin X sin - t cos ^j) dt 



where limits (-T, T) depend upon Sj^, and This integral is 
calculated using trapezoidal rule with step size (As) for 
different values of L and K. 
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method [see Eamakriahaa et.al.(19S5) for details] is as 
follows. Ue give the deaeriptioa for the system (5.3.1) only. 
With obvious notational modifications the method applies to 
(5.3.2) also. 

Let oIq be an initial approximation to the solution which 
may either be computed by some other method or guessed. Then 
this iterative method, assuming A to be of order mxn, follows 
the steps outlined below 

(i) Start with initial vector . 

0 0 

(ii) Compute the residual £q = e - A a . 

(iii) Project the residual computed in (ii) on the j-th 
column, a^ , of the matrix A and compute its component along 
the j-th column of A, i.e., 

k k ,p.k -s//" s'v 

r. = r. . - C (r. a.) / (a., a.)> a. 

— j ^ — j— “"j ”3 “"J ""j 

alter the j-th entry in the (k+l)-th approximate solution 
vector by adding the scalar component of the residual 
computed in (ii) to it, i.e., 

k+1 k,,,k N// 

otj = + c CEj_]:. §j) y 

(iv) Repeat step (iii) Pot all j=l,...,n . 

(v) Continue this process till a satisfactory solution is 


obtained . 



Using 200 iterations of this method we estimated q. (and 
. With this estimated a substituted in equation CS-3.1) the 
ESE in the reconstruction of the last 20 images of sample 1 
were predicted and compared with calculated errors. Ue found 
a maximum of 19% error in the prediction, while average error 
in prediction was 6% Table 5.1 summerizes these results. 

Similarly, with estimated^ substituted in (5.3.2). the 
predicted ESE in reconstruction of the last 20 images of 
sample 1 were obtained and compared with actual calculated 
error. The maximum error in prediction was found to be about 
28% and average error was 6.2% Table 5.2 summerizes these 


f indings . 
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TABLE 5.1 


ERROR 

PREDICTED ERROR 


2.08566300 

2.09481700 

0.00438894 

1.97690200 

1.91714200 

0.03022904 

2.04099200 

2.10144500 

0.02961927 

2.08814500 

2 .12216800 

0.01629332 

2 . 16983300 

2.07948200 

0.04163955 

2.01479800 

2.21776700 

0.10073940 

2 .41976900 

2.12497600 

0.12182710 

2.17472000 

2.08310500 

0.04212713 

2.09823000 

2.13405600 

0.01707449 

1 .96724400 

2.35344900 

0.19631780 

2.10090700 

2.05008400 

0.02419119 

2 .10284600 

2.05359300 

0.02342184 

2.01468100 

2.16347300 

0.07385401 

2.06116200 

2.34271400 

0.13659870 

2.03491500 

2.36155800 

0.16051950 

2.00582500 

2.24792500 

0.12069830 

2.03574300 

1.99760200 

0.01873579 

2.08930900 

2.09993100 

0.00508410 

2.18426200 

2.16573700 

0.00848129 

2.07897900 

2.17771100 

0.04749073 

MAXIMUM = .1963, 

MINIMUM = .0050, 

AVERAGE = .0609 
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TABLE S.B 


ERROR 

PREDICTED ERROR 

RELATIVE ERROR 

2.08566300 

2.23639500 

0.07227033 

1.97690200 

1.98593900 

0.00457148 

2.04099200 

2.15445300 

0.05559099 

2.0S814500 

2.09500600 

0.00328545 

2.16983300 

2.10870000 

0.02817424 

2.01479800 

2.15251300 

0.06835194 

2.41976900 

1.72377900 

0.28762670 

2 .17472000 

2.01027300 

0.07561736 

2.09823000 

1.98472000 

0.05409770 

1 .96724400 

2.27539400 

0.15664050 

2.10090700 

2.14389300 

0.02046053 

2.10284600 

2.03565000 

0.03195485 

2.01468100 

2.04952900 

0.01729678 

2.06116200 

2.13318000 

0.03494067 

2.03491500 

2.15644000 

0.05971985 

2.00582500 

2.28223400 

0.13780310 

2.03574300 

1.95415300 

0.04007882 

2.08930900 

2.18738000 

0.04693926 

2.18426200 

2.09159400 

0.04242534 

2.07897900 

2.11380400 

0.01675082 

MAXIMUM = .2817, 

MINIMUM = .0045, 

AVERAGE = .0627 
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S.4, AN EMPIRICAL VALIDATION OF c.i -ESTIMATE 

e 


In section 3.3 we derived an w -estimate of ESE which 

c 

depends on , namely 


C ) < 


B 


k 

k “c 




The expression implies that order of error in reconstruction 
decreases and hence the accuracy in reconstruction increases 
with the increase in the smoothness of the covariance kernel, 
provided K « *^2p' 

Now, as the following derivation shows, we note that the 
use of 'sine’ window implies p = 1. 

From equation (2.3.1) with polar co-ordinates, we have. 


NCR) 


2n CO 

„ , . -2inRr cos 4' j j ^ 
K(r)6 ^rdrd<^ 

0 0 


Differentiating it k-times and putting E = 0 , we get 


U (0) = 0, for k odd, and 


ee 


2n 


IS^ (0) = C2tT)^ C-l)^'^^ 


r^'^^ K (r) dr ( 


COS # d # > , 


for even k. 
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Thus , U^^^CO) = 0 and , 


( 0 ) 



for 'Sine’ window, it 
and hence p = 1. 

For an empirical 


K (r) dr 0 

follows that the 

validation of a 


corresponding K « K^, 


estimate we first 


c 

selected at random eight images from sample 1. For discussion 

in this section we refer to this set of S images selected 

from sample 1 as sample 9 . The reconstruction was done by CBP 

algorithm for all eight images in each one of samples 2 to 9, 

using ’Sine' window with equal to 16 and 32. For samples 5 

and 9 for each one of the total 16 images, the reconstruction 

was done for R =64 also. The number of rays N are chosen 
c • y 

as per sampling criterion and the number of views N 

& 

according to the following equation [see Natterer (1986)], 


N - 1 

y 


to reduce the discretization error in backproj action 
integral. The values of N^, total number of rays and 
total number of views for different values of R^ are given in 


Tabl e 5.3. 



TABLE 5.3 : RECONSTRUCTION PARAMETERS 


Ray-spacing, 7s 

R 

c 

Number of views 

number of Raya 

2/64 

16 

100 

64 

2/128 

32 

200 

128 

2/256 

64 

400 

256 


Not fi : 

(a) Sampling R^ = ( 1/2 7s ) 

(b) Optimal number of views = (( number of rays-l)/2 ) n 
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The errors reported are as follows: 


~ MP Z 2 ^ ^ ~ ^ I 

i j 

= -w2 2« " - - c-ij) >' 

i j 

ERMax = I i Cx. .3 - z (x^,) \ 


E = 


MP 


n I 

i j 


'CD 


Cx. .) - 
ij 


CD , , 

z (x. . ) 
ij 


where superscript CD represents the coded values of function 
in 32 gray levels. 

The error ERl represents the usual error, ERS 

represents the ESE, ERMax is the fflaximum error in a single 
pixel reconstruction, and E represents an index of the visual 
quality of the reconstructed image. For all error 
computations the region was taken to he a circle of unit 
radius. The normalised error corresponding to ERl, ERS, and E 
are given by, 

ERIN = C CERl) MP 5 / C 2 2 ^ ^ 

i j 

ERSN = C (ERS) MP ) / C 2 2 ^ 

i j 

EN = C (E) HP ) / C 2 2 ^ 

i i 
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TABLE 5. S : SUMMARY OF ERRORS IM RECONSTRUCTION 
OF RANDOM IMAGES FROM SAMPLE 3 


PIC 

R 

c 

ERl 

ERS 

ERMax 

E 

1 

16 

1.60480 

4.58385 

14.44093 

4.64804 


32 

0.54412 

0.46115 

6.87797 

0.67109 

2 

16 

1.57853 

4.50921 

mn 

4.60985 


32 

0.54594 

0.46128 


0.68876 

3 

16 

1.54895 

4.26382 

12.45509 

4.32481 


32 

0.52774 

0.42986 

6.55274 

0.64331 

4 

16 

1.57631 

4.40098 

12.85159 

4.55492 


32 

0.53904 

0.43228 

5.88706 

0.64078 

5 

16 

1.56878 

4.37426 

10.82437 

4.45170 


32 

0.53799 

0.43413 

5.68387 

0.64394 

6 

16 

1.58666 

4.47178 

13.26219 

4.53946 


32 

0.54184 

0.44864 

6.43127 

0.65309 

7 

16 

1.57031 

4.35021 

11.83567 

4.45234 


32 

0.53137 

0.42996 

6.17244 

0.64173 

8 

16 

1.54074 

4.28595 

12.16716 

4.36837 


32 

0.51318 

0.39728 

6.02370 

0.58965 
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TABLE 5.6 : SUMMARY OF ERRORS IN RECONSTRUCTION 
OF RANDOM IMAGES FROM SAMPLE 4 


PIC 

R 

c 

ERl 

ERS 

ERMax 

E 

1 

16 

1.74229 

5.25921 

14.40790 

5-36963 


32 

0.57479 

0.50102 

6.79373 

0.70739 

2 

16 

1.71181 

5.2006? 

13.59357 

5.30492 


32 

0.57437 

0.50287 

6.17243 

0-74211 

3 

16 

1 . 68099 

4.91581 

12.68245 

5.02557 


32 

0.55939 

0.46859 

6.94053 

0.68497 

4 

16 

1.71971 

5.10323 

12.94719 

5.15404 


32 

0.57247 

0.47359 

6.00125 


5 

16 

1.70936 

5.02695 

11.61611 

5.10669 


32 

0.56798 

0.47564 

6.04265 

0.69602 

6 

16 

1.72677 

5.15912 

14.12736 

5-23674 


32 

0.57384 

0.49471 

6.85965 

0.73990 

7 

16 

1.71373 

5.06646 

11.96473 

5.15530 


32 

0.56787 

0.47624 

6.15964 

0.70770 

8 

16 

1.67886 

4.94280 

12.85695 

5-04324 


32 

0.54414 

0.43551 

5.98677 

0. 64994 
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TABLE 5.7 : SUMMARY OF ERRORS IK RECOKSTRUCTION 
OF RANDOM IMAGES FROM SAMPLE 5 


PIC 

R 

c 

ERl 

ERS 

ERMax 

E 


16 

1.69169 

4.81965 

11.11550 

4.92456 

1 

32 

0.53339 

0.43641 

5.63172 

0.66162 


64 

0.31896 

0.12656 

1.48591 

0.21717 


16 

1-70300 

4.95327 

12.51725 

5-06092 

2 

32 

0.5420? 

0.44455 

5.66347 

0.67361 


64 

0.32468 

0.13217 

1.73761 

0.22348 


16 

1.73848 

4.92837 

11.37895 

4.99053 

3 

32 

0.54476 

0.43898 

6-52774 

0-64583 


64 

0.32535 

0.13312 

1.88581 

0.21622 


16 

1.71720 

4.88942 

11.56064 

4-98864 

4 

32 

0.53903 

0.42903 

5.52768 

0.62689 


64 

0.32080 

0-12840 

1.65534 

0.21559 


16 

1.71076 

4.91810 

10.47233 

5.04924 

5 

32 

0.53676 

0.43074 

5.71330 

0.64047 


64 

0.31552 

0.12637 

1.61744 

0-20960 
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TABLE 5 . 7 (Cofitd. . . . ) 


PIC 

R 

c 

ERl 

ERS 

ERHax 

1 


16 

■ 

4.81116 

10.96620 

4.88573 

6 

32 

0.5314S 

0.42293 

5.60305 

0.63068 


64 

0.32138 

0.12679 

1.72970 

0.21465 


16 

1.71673 

4.90204 

11.31615 


7 

32 

0.53072 

0.43050 

6.01982 



64 

0.31942 

0.12699 

1.67328 

0.20581 
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TABLE 5.S : SUMMARY OF ERRORS IN RECONSTRUCTION 




OF RANDOM 

IMAGES 

FROM SAMPLE 6 


PIC 

R 

c 

ERl 

ERS 

ERMax 

E 

1 

16 

2.29068 

8.13229 

13.33023 

8.15972 


32 

0.67101 

0.64205 

6.67258 

0 . 90720 

2 

16 

2.28511 

8.24757 

14.09712 

8.35006 


32 

0.68502 

0.65850 

6.17501 

0.94097 

3 

16 

2.33983 

8.34802 

13.72709 

8.40278 


32 

0.68784 

0.65192 

7-67429 

0.92077 

4 

16 

2.32022 

8.30365 

13.95919 

8.37500 


32 

0.68889 

0.64639 

6.21361 

0-90246 

5 

16 

2.29941 

8.28919 

12.49980 

8.44634 


32 

0.68162 

0.64259 

6.70382 

0.92393 

6 

16 

2.28556 

8.16798 

12.63069 

8.29640 


32 

0.67729 

0.63584 

6.51307 

0.91477 

7 

16 

2.29812 

8.26638 

14.21908 

8-28725 


32 

0.67548 

0.64520 

6.64463 

0.93813 

8 

16 

2.28029 

8.01864 

12.16393 

8.13920 


32 

0.66887. 

0.60333 

6.33409 

0.87689 
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TABLE S. S ; SUMMARY OF ERRORS IN RECONSTRUCTION 
OF RANDOM IMAGES FROM SAMPLE 7 


PIC 

R 

c 

ERl 

ERS 

ERMax 

E 

1 

16 

2.72255 

11.18504 

14.51111 

11.28946 


32 

0.77759 

0.82865 

7.36782 

1.15877 

2 

16 

2.70707 

11.29198 

15.64118 

11.38352 


32 

0.79302 

0.85694 

6.86103 

1.20044 

3 

16 

2.77559 

11.50256 

14.66276 

11.58807 


32 

0.79561 

0.84559 

8.58786 

1.19413 

4 

16 

2.75426 

11.43549 

15.51450 

11.56755 


32 

0.80029 

0.84789 

6.7007? 

1.18624 

5 

16 

2.73954 

11.45025 

13.39005 

11.49811 


32 

0.7S9S3 

0.83917 

7.20628 

1.15783 

6 

16 

2 .72089 

11.30927 

14.42341 

11.42929 


32 

0.78656 

0.83274 

7.40941 

1.15657 

7 

16 

2.73190 

11.38676 

15.40907 

11.48927 


32 

0.7S7S8 

0.84354 

7.51570 

1.16477 

8 

16 

2.69854 

11.02025 

12.76122 

11.10006 


32 

0.77123 

0.78056 

7.05621 

1.08744 
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TABLE S. lO ; SUMMARY OF ERRORS IN RECONSTRUCTION 
OF RANDOM IMAGES FROM SAMPLE 8 


PIC 

R 

c 

ERl 

ERS 

ERMax 

1 

1 

16 

2.97979 

13.26143 

15.35691 

13 . 31660 


32 

0.84243 

0.95812 

7.71013 

1.33112 

2 

16 

2.95555 

13.32545 

16.28948 

13.38794 


32 

0.85628 

0.98611 

7.27841 

1.35038 

3 

16 

3.03090 

13.63298 

15.30843 

13.76894 


32 

0.86083 

0.97917 

9.01085 

1.37374 

4 

16 

3.01287 

13.58625 

16.63249 

13.66635 


32 

0.86702 

0.98326 

7.11649 

1.36774 

5 

16 

2-99357 

13.54150 

14.17572 

13.56503 


32 

0.85329 

0-96535 

7.67799 

1.32891 

6 

16 

2.97468 

13.43073 

15.36428 

13.42550 


32 

0.85440 

0.96795 

7.87311 

1.32860 

7 

16 

2.98069 

13.46921 

16.27425 

13.51894 


32 

0.86426 

0.97622 

7.69080 

1.32544 

8 

16 

2.96074 

13.04525 

13.64192 

13.10922 


32 

0.83340 

0.90135 

7.19785 

1.23674 
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TABLE S.ll ; SUMMARY OF ERRORS IN RECONSTRUCTION 
OF RANDOM IMAGES FROM SAMPLE 9 


PIC 

R 

c 

ERl 

ERS 

ERMax 

E 


16 

0-91947 

2.17561 

9.61702 

2.24116 

1 

32 

0.389S4 

0.27742 

4.84390 

0.43529 


64 

0.22451 

0.06636 

1.47740 

0.07670 


16 

0.91210 

2.18695 

9.80218 

2.29072 

2 

32 

0.39307 

0.28037 

4.83003 

0.43561 


64 

0.22779 

0.06808 

1.48973 

0.07860 


16 

0.90719 

2.12095 

10.47944 

2.18119 

3 

32 

0.39102 

0.27301 

5.48529 

0.41667 


64 

0.22703 

0.06761 

1.35773 

0.07860 


16 

0.91799 

2.14057 

9.38185 

2.22001 

4 

32 

Q.3S794 

0.26748 

5.02848 

0.40972 


64 

0.22472 

0.06597 

1.38346 

0.07355 


16 

0.91086 

2-13622 

8.63433 

2.23264 

5 

32 

0.38858 

0.26859 

4.78136 

0.41383 


64 

0.22516 

0.06598 

1.30920 

0.07607 
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TABLE 5.13 : SUMMARY OF NORMALIZED ERRORS IN RECONSTRUCTION 
OF RANDOM IMAGES FROM SAMPLE 3 


PIC 

R 

c 

■ 

ERSN 

EN 

1 

16 

0.10476 

0.13507 

0.30342 


32 

0.03552 

0.04871 


2 

16 

0.10198 

0.13258 

0.29782 


32 

0.03527 

0.04822 

0.04450 

3 

16 

0.09850 

0.12736 

0.27502 


32 

0.03356 

0.04598 

0.04091 

4 

16 

0.10235 

0.13192 

0.29575 


32 

0.03500 

0.04701 

0.04161 

5 

16 

0.10136 

0.13072 

0.28763 


32 

0.03476 

0.04683 

0.04161 

6 

16 

0.10287 

0.13227 

0.29431 


32 

0.03513 

0.04764 

0.04234 

I 

16 

0.10222 

0.13131 

0.28983 


32 

0.03459 

0.04694 

0.04177 

8 

16 

0.10220 

0.13288 

0.28976 


32 

0.03404 

0.04600 

0.03911 
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table 5.14, ; SUMHARY OF NORMALIZED ERRORS IN RECONSTRUCTION 
OF RANDOM IMAGES FROM SAMPLE 4 


FIC 

R 

c 

ERIN 

ERSN 

EN 

1 

16 

0.11370 

0.14405 

0.35042 


32 

0.03751 

0.05056 

0.04616 

2 

16 

0.11057 

0.14180 

0.34266 


32 

0.03710 

0.05014 

0.04793 

3 

16 

0.10689 

0.13624 

0.31956 


32 

0.03557 

0.04783 

0.04356 

4 

16 

0.11172 

0.14158 

0.33483. 


32 

0.03719 

0.04904 

0.04501 

5 

16 

0.11045 

0.13958 

0.32997 


32 

0.03670 

0.04882 

0.04497 

6 

16 

0.11200 

0.14153 

0.33966 


32 

0.03722 

0.04983 

0.04799 

7 

16 

0.11157 

0.14111 

0.33563 


32 

0.03697 

0.04919 

0.04607 

S 

16 

0.11135 

0.14207 

0.33449 


32 

0.03609 

0.04795 

0.04311 
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TABLE S. 15 : SUMMARY OF NORMALIZED ERRORS IN RECONSTRUCTION 
OF RANDOM IMAGES FROM SAMPLE 5 


PIC 

R 

c 

ERIN 

ERSN 

EN 


16 

0-11053 

0.14059 

0.32176 

1 

32 

0.03435 

0.04310 

0.04323 


64 

0.02034 

0.02590 

0.01419 


16 

0.10999 

0.14089 

0.32636 

2 

32 

0.03501 

0.04799 

0.04351 


64 

0.02097 

0.02617 



16 

0.11061 

0.13863 

0.31752 

3 

32 

0.03466 

0.04705 

0.04109 


64 

0.02070 

0.02591 

0.01376 


16 

0.11150 

0.14034 

0.32392 

4 

32 

0.03500 

0.04744 

0.04070 


64 

0.02033 

0.02595 

0.01400 


16 

0.11050 

0.14042 

0.32614 

5 

32 

0.03467 

0.04725 

0.04137 


64 

0.02038 

0.02559 

0.01354 




RIW 


ERSN 


EH 


6 

16 

32 

64 

0.10998 

0.03448 

0.02085 

0.13944 

0.04701 

0-02574 

0.31696 

0.04092 

0.01393 


16 

0.11163 

0-14107 

0.32422 

7 

32 

0.03451 

0.04754 

0.04206 


64 

0.02077 

0.02582 

0.01338 


16 

0.11269 

0.14188 

0.32263 

S 

32 

0.03508 

0.04705 

0.04101 


64 


0.02061 


0-02563 


0-01284 
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TABLE 5.16 : SUMMARY OF NORMALIZED ERRORS IN RECONSTRUCTION 
OF RANDOM IMAGES FROM SAMPLE 6 


PIC 

R 

c 

ERIN 

ERSN 

In 

1 

16 

0.14966 

0.17917 

0.53311 


32 

0.04384 

0.05724 

0.05927 

2 

16 

0.14761 

0.17842 

0.5393S 


32 

0.04425 

0.05733 

0.06078 

3 

16 

0.14879 

0.17709 

1 


32 

0.04374 

0.05627 

0.05855 

4 

16 

0.15072 

0.18022 

0.54403 


32 

0.04475 

0.05717 

0.05862 

5 

16 

0.14860 

0.17897 

0.54585 


32 

0.04405 

0.05666 

0.05971 

6 

16 

0.14821 

0.17804 

0.53799 


32 

0.04392 

0.05648 

0.05932 

7 

16 

0.14956 

0.17989 

0.53933 


32 

0.04396 

0.05715 

0.06105 

a 

16 

0.15123 

0.18072 

0.53980 


32 

0.04436 

0.05637 

0.05816 
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TABLE S. 17 : SUMMARY OF NORMALIZED ERRORS IN RECONSTRUCTION 
OF RANDOM IMAGES FROM SAMPLE 7 


PIC 

R 

c 

ERIN 

ERSN 

In 

1 

16 

0.17762 

0.20630 

0.73653 


32 

0.05073 

0.06385 

0.07560 

2 

16 

0.17471 

0.20510 

0.73467 


32 

0.05118 

0.06425 

0.07747 

3 

16 

0.17649 

0.20454 

0.73684 


32 

0.05059 

0.06306 

0.07593 

4 

16 

0. 17S86 

0.20790 

0.75119 


32 

0.05197 

0.06437 

0.07703 

5 

16 

0.17686 

0.20662 

0.74230 


32 

0.05099 

0.06360 

0.07475 

6 

16 

0.17649 

0.20589 



32 

0.05102 

0.06353 

0.07502 

7 

16 

0.17767 

0.20741 

0.74721 


32 

0.05124 

0.06419 

0.07575 

a 

16 

0.17901 

0.20829 

0.73633 


32 

0.05116 

0.06303 

0.07214 




134 


TABLE 5.18 : SUMMARY OF NORMALIZED ERRORS IN RECONSTRUCTION 
OF RANDOM IMAGES FROM SAMPLE 8 


PIC 

R 

c 

ERIN 

ERSN 

EN 


1 

16 

0.19433 

0.22208 

0.86846 



32 

0.05494 

0. 06788 

0.08681 


2 

16 

0.19077 

0.22039 

0.86414 



32 

0.05527 

0.06817 

0.08716 


3 

16 

0.19291 

0.22050 

0.87636 



32 

0.05479 

0.06719 

0.08744 


4 

16 

0.19571 

0-22414 

0.88774 



32 

0.05632 

0-06856 

0.08885 


5 

16 

0.19348 

0.22244 

0.87673 



32 

0.05515 

0-06753 

0.08589 


6 

16 

0.19281 

0.22168 

0.87020 



32 

0.05538 

0.06767 

0.08612 


7 

16 

0.19393 

0.22316 

0.87957 



32 

0.05558 

0.06831 

0.08624 


8 

16 

0.19583 

0.22414 

0.87001 



32 

0.05531 

0.06699 

0.08208 
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table S. 19 : SUMMARY OF NORMALIZED ERRORS IN RECONSTRUCTION 
OF RANDOM IMAGES FROM SAMPLE 9 


PIC 

R 

c 

ERIN 

ERSN 

EN 


16 

0.06012 

0.09553 

0.14654 

1 

32 

0.02549 

0.03879 



64 

0.0146S 

0.01897 

0.00502 


16 

0.05S94 

0.09464 

0.14803 

2 

32 

0.02540 

0.03853 

0.02815 


64 

0.01472 

0.01899 

||n^^^H|||||||| 


16 

0.05770 

0.09187 

0.13873 

3 

32 

0.02487 

0.03748 

0.02650 


64 

0.01444 

0.01865 

0.00500 


16 

0.05956 

0.09414 

0.14404 

4 

32 

0.02517 

0.03784 

0.02658 


64 

0.01458 

0.01879 

0.00477 


16 

0.05886 

0.09358 

0.14427 

5 

32 

0.02511 

0.03773 

0.02674 


64 

0.01455 

0.01870 

0.00492 
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Fig. 5.3 


Reconstruction of Fig 


5.1 for R 

c 


Fig. 5.4 



Reconstruction of Fig. 5.2 for 


= 16. 


16 
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Fig. 5. 


Fig. 5. 




• :-r^l;r>.5i Knl 'jCtU ^~■ 

Reconstruction of Fig. 5.1 for R 

c 


32 


• ^ ‘I** 


S Reconstruction of Fig. 5.2 for R^ - 32. 
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Fig. 5,7 


Reconstruction .of Fig. 


5.1 for R 

c 



= 64. 


Fig. 5.8 


64 
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CHAPTER 6 

THE OPTIMAL FILTER 


6.1 INTRODUCTION 


In chapters 1 and 2 we have discussed the CBP Hiethod and 
it's discrete implementation. This method idealizes the 
projection data as line integrals, but in actual practice the 
projection data are integrals over strips o£ finite width. In 
the case of Positron Em^isssion Tomography CPET> , the data 
collection procedure is entirely different. However, the PET 
problem gets reduced to the CT problem, using transmission 
tomography, after a suitable attenuation compensation has 
been made [see, e.g., Vardi, Shepp and Kaufman (1985)1 . 

Thus, as shown in section 2.2, in implementing CBP 
algorithm, approximations are introduced both while taking 
care of data collection procedure, as well as while 
performing discrete backpro jection operation. In the CBP 
method one would wish for a convolving function (also called 
filter in this chapter) q which does the data-noise reduction 
and compensates for the discretized model, as well as the 
interpolations made. The question, however, 


is whether a 
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single convolving function could be chosen to do all these 
tricks simult enously . 

The present practice seems to be that of experimenting 
with various analytical convolving functions [see 
Ramachandran and Lakshminarayanan (1971) and Shepp and Logan 
(1974)] on a set of standard phantoms and finally using the 
one which works beat for the type of crosa“sect ion under 
study. 

In this chapter we propose a method of derivation of a 
convolving function which optimally takes care of (i) the 
stochastic process giving rise to object image zCx,y), (ii) 
the data collection geometry, (iii) the cross-sectional 
discretization scheme, and (iv) the data noise, while at the 
same time minimizing a given expected weighted sum of squared 
error in pixel reconstruction. 

6.S THEORETICAL FQRMULATIOM 

Ue assume that the discretization of object, whose 
intensity or cross-sectional distribution at point (x,y) is 
given by z(x,y), is done on an NxN grid of pixels. Ue also 
assume that object distribution is limited to a unit circle 
with origin as centre. 

z = Cz Cx .y), i,i=i. .N > 


Let 
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denote the tf-vector (written in matrix form) of discretized 
object distribution. Let the KxL matrix P = C P 
denote the matrix of projection data where K and L are the 
total numbers of rays and views, respectively, for which data 
is available and is "the projection along k~th ray 

and 1 -th view . 

For the reconstruction z Cx. , x 3 of z fx , x > by 
CBP method, equations (2.2.1 - 2.2.3) can be rewritten as 
follows : 

L 

z Cx . ,x^ .> = (As) p Cx cos & + X sin e , e > 

I =1 

L K 

= (A6) 2 CAs’) 2 (As) 2 <1 t®k'" ®k>- 

1=1 k' k=l 

I Cx . cos &. + x„ .sin 6 - s ,> ...(6.2.1) 

'•it i 24 Ik 

L K 

= (A©)(As)(As-) 2 2 2 * 

l =1 k' k =1 

I Cx .cos a + X sin e - s ,> 

L K 

= (Aa)(As)(As') 2 2 Z ' ®k' 

I =1 k' k =1 


I Cx^j^cos + x^^sin a^^- Sj^,> 
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L K 

' J J (Aa-) I 3^, 

1=1 k=l k' 

I Cx^.cos - Sj^,).q Cs^^) 

K L 

= (Ae)CAs) 2 2 ®l" ®k’ ^ <®k>- 

k =1 I =1 

... ( 6 . 2 . 2 ) 

where denote the modified projection data. 

Thus in matrix form (6.2.2) can be written as 

=z=Xq (6.2.3) 

where z is N^-vector of reconstructed object distribution, X 

is an N^x K matrix and the K-vector q is the convolving 

funct ion (filter). As the k-th element <l(Sj^) of q is taken 

to be independent of the particular view and depends on the 

magnitude of alone, it will be denoted by and will be 

called k-th filter coefficient. 

The error in (i,j)-th pixel reconstruction is 

z(x ,y ) - z(x ,y.) and the weighted total expected 
>■ j '■ i 

reconstruction error is 

» = [ Cz “ Cz - 2)] ....,(6.2.4) 

where T is a positive semi-definite matrix representing , the 
.relative weights given to the pixels. In. fiainifaizing the 
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expected error 1 we get the normal equations 

^ ® (6.2.5) 

where Q. C5£ T X) and R = (X'T z) with X-matrix as in 

equation (6.2.3). As is well known (6.2.5) is a consistent 
system, and any solution of the system (6.2.5) gives an 
optimal filter. In cases when Q is deficient the q of 
interest is the Moore-Penroae solution of (6.2.5). Hote that 
Q is positive semi-definite and would be positive definite in 
most cases of practical interest. 

When the parameters of the process ^ are unknown we 

propose the following method to estimate the optimal filter 
from the ensemble. 

6.3 ESTIMATION OF OPTIMAL FILTER 


Let z^,...,z^ be a random sample (consisting of 
phantoms, in practice ) from X and let X^,,.,,X^ be the 
corresponding X~matrices determined from projection data on 
the machine. From the law of large numbers, since z^and X^ , 
i= l,...,n are independently and identically distributed, it 
follows that 


n ^ 

O = - y X.'T X. and R = ^ T ^ 'T z 
nnZi- '- nn/ii- ‘• 


(6.3.1) 


i =1 


i =1 


are strongly consistent estimators of G and R, 


respectively. 
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+ + 

= ft E i. d R = a w o 1 

n n n ^ ^opt w.p.i 

i. R w.p.i, and d is unique it follows that the 

convergence w.p.i is along the whole sequence rather than 

just along a subsequence. In such a situation where a choice 

of r is not readily estimable, one may use a sufficiently 

small « and truncate those xf’^’s which are smaller than e 

1 - 

and estimate them as 0 . 

(II) APPROXIMATE EXPECTATION HIHIHIZATIQN (AEM) : Let 

z , . . . ,z be reconstructions of z z . In AEM we 

si r» in 

estimate q by minimizing the error 

Opt 


e 

r» 


i \ [C z - ^ )'T (z - z>] 

Tt j£j It. 1.1. 


... (6.3.3) 


Since by law of large numbers, e y * w.p.i, for 

n 

suf f ieciently large n a q which minimizes e , will 
approximately minimize the errror S in (6.2.5). 


6.4 NUMERICAL IMPLEMENTATION AND CQMPARIS<»IS 


For the purpose of numerical implemetat ion of the 
optimal filter, the sample taken is of 17 40x40 digitized 
pictures from Gonzale*(1977), Munshi(19SS)and Eaatogi (1988). 
Projection data for these object distributions (pictures) 


were obtained using: 
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(a) Transiftisaion (CAT) method and line integrals, 

(b) 1 ransmi as ion (CAT) method and strip integrals, 

(c) Emission (PET) method and line integrals, 
and 

(d) Emission (PET) method and strip integrals, 


thus giving 

rise to four 

modes of 

data collection. 

The 

data 

acquisition 

was carried 

out for 40 

views and 40 

rays in 

the 

case of CAT 

and with 40 

detectors 

in PET case. 

As 

shown in 

Figure 6.1, 

calculations 

of strip integrals 

in 

PET 

ease 

involve strips of width 


cos(<9+^). where, 

Z6 is 

the 


angle the two detectors defining the atrip make at the centre 
and 14 > is the angle the two adjoining detectors make at the 
centre. It is assumed that the detectors are fixed at equal 
distances around the object on a circle in the same plane. 
All calculations for obtaining the optimal filter were 
carried out with matrix T as identity and Gaussian 
elimination method was used for solving the resulting system 
of linear equations. The comparisons have been done with 
standard filters. [Ramachandran and Lakshminarayanan (1971), 
Shepp and Logan (1974), and CSI filter of Rathore et.al. 
(1988), RKS filter of Rathore( 1987 ) ] with respect to the 
error criterion 

EL = [ J ( i. / HP 

2, ^ L 

2 and z are the actual object function and the 
ij ^ M 
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recoristt-ucted object ftiftetion values in Ci,j)-th pixel and MP 
ia the total nuiftber of pixels used in the suianiatiori. 

Results (values of in the noise free case are given 

in Tables 6.1~6.4. These results for test case pictures are 
shown in Figures 6.2-6.10, for CAT line integral, i.e., case 
(a), Figures 6.11-6.19 for case (b). Figures 6. 20-6. 28 for 
case (c) and Figures 6.29-6.37 for case (d). 

In order to simulate a noisy case, independent and 

identically distributed pseudo random variables having a 

uniform distribution over the interval (-7 . 5 , 7.5) were 

generated, projection data for noise alone were obtained and 

the resulting noise was then added to the noise-free 

projection data, and from this noisy, data reconstructions 

were made. The EL -errors in this case are given in Tables 

2 

6. 5-6. 8. The corresponding results for test case pictures are 
sown in Figures 6.38-6.46 for case (a). Figures 6.47-6.55 for 
case (b) , Figures 6.56-6.64 for ease (c) and Figures 
6.65-6.73 for ease (d). 


6.5 CONCLUSIONS 

From the tables of EL it is clear that the proposed 
filter cofflpares very faLVoraibly with afiy o£ the other filters 

used. For certain pictures, in fact, the reduction in error is 
upto 20 percent, in noise-free case and upto 30 percent, in 
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tioiey case. 

lii£: followltig features of the proposed, filter are worth 
noting: 

(1) Since the class of pictures on which it was tested 
is quite heterogeneous, the derivation gives a filter that is 
quite robust, as the results show. 

(2) A comparative study of noisy and noise-free eases 
shows that as compared with the standard filters, the optimal 
filter in noisy environment is indeed effective in noise 
reduction. Thus in the noisy ease, without prefiltering the 
projection data, the performance of CBP with this optimal 
filter is expected to be better than that of other methods 
like ART. 

(3) Provisions for further improvement of the filter are 
also present. As new cross-sectional distributions (phantoms) 
become available, their information can also be incorporated 
in deriving the improved filter. 




FIG^.l Showir^ width of strip in PET strip. 
(D2, D18) is detector strip. 

C is the centre. 


1S3 


table 6.1: ERROR IN CAT LINE INTEGRAL (NOISE-FREE) CASE 


PIC/FILTER OPT 

RAM 

RKS 

SHE 

CSI 

LIN 

3.406 

3.696 

3.870 

3.657 

3.754 

FAC 

3.304 

4.499 

4.695 

4.412 

4.521 

BRA 

4.661 

4.732 

4.964 

4.661 

4.697 

STA 

3.146 

3.182 

3.326 

3.146 

3.206 

COR 

13.672 

14.016 

14.634 

13.809 

13.878 

SAT 

3.818 

5.325 

5.551 

5.250 

5.401 

THO 

4.077 

4.174 

4.402 

4.093 

4.142 

SMO 

0.941 

1.041 

1.083 

1.026 

1.026 

PLA 

2.913 

3.141 

3.299 

3.071 

3.089 

BOY 

3.574 

3.437 

3.613 

3.398 

3.496 

MON 

2.344 

3.177 

3.330 

3.133 

3.221 

UOR 

3.135 

3.161 

3.315 

3.119 

3.151 

CHA 

5.229 

6.270 

6.537 

6.190 

6.350 

FIN 

2.823 

3.751 

3.872 

3.710 

3.791 

CRO 

3.445 

3.531 

3.704 

3.488 

3.531 

JET 

4.067 

4.235 

4.423 

4.179 

4.217 

GEO 

2.372 

3.276 

3.374 

3.227 

3.349 


AVERAGE 3.937 4.391 4.5S8 4. 328 4.401 
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TABLE 6.2: ERROR IN CATSTRIP INTEGRAL (NOISE FREE) CASE 


PIC/FILTER 

OPT 

RAM 

RKS 

SHE 

CSI 

LIN 

3.464 

3.793 

3.619 

3.657 

3.735 

FAC 

3.347 

4.586 

4.390 

4.455 

4.521 

BRA 

4.714 

4.910 

4.643 

4.697 

4.697 

STA 

3.170 

3.254 

3.122 

3.146 

3.194 

COR : 

L3.649 

14.313 

13.695 

13.809 

13.787 

SAT 

3.S6S 

5.426 

5.225 

5.250 

5.401 

THO 

4.126 

4.321 

4.061 

4.126 

4.142 

SMO 

0.941 

1.069 

1.012 

1.026 

1 .026 

PLA 

2.966 

3.247 

3.054 

3.106 

3.089 

BOY 

3.633 

3.535 

3.359 

3.398 

3.476 

MON 

2.410 

3.243 

3.111 

3.133 

3-221 

UOR 

3.135 

3.246 

3.087 

3.113 

3-135 

CHA 

5.256 

6.377 

6.137 

6.190 

6.323 

FIN 

2.843 

3.771 

3.690 

3.690 

3.791 

CRO 

3.463 

3.643 

3.463 

3.497 

3.514 

JET 

4.104 

4.348 

4.160 

4.198 

4.217 

GEO 

2.494 

3.301 

3.203 

3.227 

3.349 

AVERAGE 

3.975 

4.493 

4.296 

4.336 

4.389 
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TABLE 6.3: L^ERROR IN PET LINE INTEGRAL (NOISE-FREE) CASE 


PIC/FILTER 

OPT 

RAM 

RKS 

SHE 

CSI 

LIN 

3.735 

4.044 

4.141 

4.122 

4.315 

FAC 

3.673 

4.912 

4.868 

4.890 

5.021 

BRA 

4.839 

4.910 

5.052 

5.017 

5.213 

STA 

3.542 

3.626 

3.698 

3.674 

3.807 

COR ; 

L5.000 

15.023 

15.092 

15.046 

15.115 

SAT 

4.346 

5.853 

5.878 

5.876 

6.104 

THO 

4.272 

4.434 

4.597 

4.548 

4.840 

SMO 

0.969 

1.112 

1.055 

1.069 

1.069 

PLA 

3.212 

3.440 

3.422 

3.422 

3.510 

BOY 

4.043 

3.808 

3.945 

3.906 

4.160 

MOM 

2.761 

3.528 

3.528 

3.506 

3.637 

UOR 

3.570 

3.628 

3.697 

3.676 

3.761 

CHA 

6.217 

7.150 

7.257 

7.230 

7.471 

FIN 

3.226 

4.194 

4.094 

4.114 

4.114 

CRO 

3.850 

3.962 

4.057 

4.031 

4.177 

JET 

4.591 

4.798 

4.798 

4.798 

4.910 

GEO 

2.958 

3.569 

3.545 

3.545 

3.667 

AVERAGE 

4.400 

4.823 

4.866 

4.851 

4.994 
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TABLE 6.6: ERROR IN PET LINE INTEGRAL (NOISY) CASE 


PIC/ FILTER 

OPT 

RAH 

RKS 

SHE 

CSI 

LIN 

4.044 

4.277 

4.238 

4.296 

4.470 

FAC 

4.10S 

5.216 

5.260 

5.194 

5.303 

BRA 

5.106 

5-159 

5.106 

5.195 

5.337 

STA 

3.963 

3.939 

3.939 

3.939 

4.011 

COR 

15.023 

15.092 

15.046 

15.115 

15.161 

SAT 

4.622 

6.054 

6.054 

6.054 

6.280 

THO 

4.5S0 

4.727 

4.629 

4.759 

4.987 

SMO 

1.768 

1.668 

1.739 

1.639 

1.597 

PLA 

3.475 

3.633 

3.668 

3.615 

3.703 

BOY 

4.160 

4.023 

3.926 

4.062 

4.258 

MON 

3.133 

3.834 

3.856 

3.834 

3.944 

UOR 

3.867 

3.877 

3.856 

3 . 893 

3.930 

CHA 

6.323 

7.230 

7.150 

7.257 

7-471 

FIN 

3.489 

4.295 

4.376 

4.255 

4.295 

CRO 

4.160 

4.195 

4.169 

4.203 

4.281 

JET 

4.798 

4.929 

4.948 

4.910 

5.004 

GEO 

3.276 

3.692 

3.765 

3-667 

3-741 


AVERAGE 4.700 5-049 5-043 5.052 5.163 





TABLE 6.7; ERROR IN CAT STRIP INTEGRAL (NOISY) CASE 


PIC/FILTER OPT RAM RKS SHE 


LIN 

4.006 

FAC 

4. IDS 

BRA 

5.213 

STA 

3.951 

COR 

13.878 

SAT 

4.421 

THO 

4.743 

SMO 

2.495 

PLA 

3.808 

BOY 

4.375 

MON 

3.396 

UOR 

3.782 

CHA 

5.576 

FIN 

3.569 

CRO 

4.177 

JET 

4.648 

GEO 

3.325 


4.509 

4.954 

5.368 

5.781 

5.390 

5.889 

4.215 

4.695 

14.176 

14.840 

6.003 

6.405 

5.019 

5.490 

3.022 

3.621 

4.194 

4.721 

4.551 

5.078 

4.229 

4.733 

4.074 

4.620 

6.697 

7.177 

4.578 

5.001 

4.479 

4.996 

4.985 

5.454 

4.254 

4.719 


4.296 

4.219 

5.151 

5.129 

5.177 

5.052 

3.987 

3.831 

13.970 

13.993 

5.828 

5.853 

4.792 

4.678 

2.666 

2-310 

3.931 

3.756 

4.277 

4.140 

3.966 

3.834 

3.808 

3.628 

6.510 

6.563 

4.376 

4.295 

4.238 

4.074 

4.760 

4.629 

4-010 

3.936 


AVERAGE 4.675 


5.279 5.775 5.044 4.936 
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TABLE 6.8: ERROR IN PET STRIP INTEGRAL (NOISY) CASE 


PIC/FILTER 

OPT 

RAM 

RKS 

SHE 

CSI 

LIN 

4.277 

4.199 

4.102 

4.238 

4.412 

FAC 

4.021 

4.521 

4.390 

4.564 

4.825 

BRA 

5.390 

5.177 

5.106 

5.230 

5.373 

STA 

3.843 

3.735 

3.686 

3.759 

3.855 

COR 15.023 

14.977 

14.932 

15.023 

15.069 

SAT 

4.773 

5.627 

5.426 

5.727 

6.079 

THO 

4.580 

4.678 

4.564 

4.710 

4.954 

SMO 

1.967 

1.639 

1.697 

1.625 

1.597 

PLA 

3.475 

3.387 

3.352 

3.387 

3.510 

BOY 

4.316 

4.043 

3.926 

4.082 

4.297 

MON 

3.374 

3.550 

3.484 

3.571 

3.747 

UOR 

4.020 

3.925 

3.914 

3.936 

3.957 

CHA 

6.670 

7.097 

6.937 

7.177 

7.444 

FIN 

3.448 

3.791 

3.731 

3.831 

3.952 

CRO 

4.264 

4.169 

4.117 

4.195 

4.298 

JET 

4.798 

4.704 

4.610 

4.741 

4.891 

GEO 

3.447 

3.569 

3.472 

3.618 

3.790 

AVERAGE 

4.805 

4.870 

4.791 

4.907 

5.062 





Reconstruction of Lincoln in CAT line integral noise 
free case (a) CSI filter, (b) RAM filter, (c) RKS 
filter, (d) SHEPP filter, (e) original picture, 

(f) optical filter. 
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Fig. 6.10 Reconstruction of Chromosomes in CAT line integral 

noiss-'f rS'B cslsb .(a) ■ CSX .filtaX/ (b) ■ RAM , rxlt©x^ (c) RKS 
filter, (d) Shepp filter, (e) original picture, 

(f ) optimal filter. 
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Pig. 6.19 Reconstruction of Chroroosomes in CAT strip integral 
noise-free case (a) CSI filter^ (b) RAI4 filter, 

(c) RKS filter, (d) Shepp filter, (e) original 
picture, (f) optimal filter. 




nstruction of Lincoln in PET 
case (a) CSI filter, (b) RW- 
.er, (d) Shepp filter, (e) oric 
ootiinal filter. , 





Reconstruction of . 
noise-free case (a 
(c) RI<S filter, (d 
(f) optiinal filter 








Reconstruction of Saturn xn PE' 
free case (a) CSI filter, lb) . 
filter, (d) Shepp filter, (e) 
(f) ootimal filter. 
















mm 






Reconstruction of Brain in 
free case (a) CSI filter, 
filter, (d) Shepp filter, 
(f) opt imal filter. 


C/l 








Reconstruction o 
noise-free case 
(c) RXS filter, 
picture, (f) opt 




Reconstruction or i- 
noise-free case (a) 
(c) RiCS filter, (d) 
picture, (f) op time 
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Reconstruction of Brain in CAT line inte 
case (a) CSI filter, (b)'RAii filter, (c) 
(d) Ghepp filter, (e) original picture, 
filter. 


Reconstruction of Statue of Liberty in CAT line 
integral noisy case (a) CSI filter, (b) RAi''. fil 
(c)RKS filter, (d) Shepp filter, (e) original 
picture, (f) optimal filter. 







■onalfisa , in CAT line 

, (b) RAi-* filter, i 

e) original picture 


















Reconstruction or j 
case (a) CSI filter 
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Reconstruction of 

noisy case CSI 
filter, (d) Shepp 


Ghromosomes 
filter, ir 
filter, (e) 


( f } ODtirnal f * 


in CAT strip integ 
RAI'-. filter, ic] R 
original picture. 







■hepp filter 





Reconstruction of Saturn in PET line inte 
case (a) CSI filter, (b) RAr'i filter, (c) 
(d) Shepp filter, (e) original picture, ( 
filter. 


Ml 














Rec on s true tion ox L rncoxn 
noisy case (a) CSI filter 
filter, (^l) Shepp filter, 
(f) optimal filter. 
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Reconstr’action of Thorax in PET strip integral noisy 
case (a) CSI filter, (b) RAM filter, {c} RKS filter, 
(d) Shepp filter, ' (e) originel picture, (f) optimal 
filter. 




Reconstruction of MonaLisa in PET strip integral 
iioisy case (a) CSi filter, (b) RAi'; filter, (c) RKS 
filter, (d) Shepp filter, (e) original picture, 

(f) optimal filter* 







Reconstruction of Chromosomes in PEC strip inte 
noisy case (a) GSI filter, (b) RAi-i filter, (c) 
filter,’ (d) Shepp filter, (e) original picture, 
(f ) optimal, filter. 
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1 

Zfi 




Cyi+CUi-Vi)/A, y^^(u^-v^)/A ), 


0 [R^£R^ 


K (u) K (v) du dv . 

Since K (x) « [0, oo) , Fubiini’s theorem ia applicable and 

thus , 


_ Cy) = 


"iS: 


A,v. ^ 


^^.r, ^ Cy + Cy - KCy) du dy 




and 




E* E^ 


Hence , 

I “r! ^ ‘^SS W 1 ^ 


^2 „2 
E E 


I ^ ^ ^ ^ ~ ^ I - 


I K (u) I I K Cy) I du dy 


^ 1 “J! J '^3 llo « 


K Cu) [ du 


E 


implying that. 


i ^ llo ~ ”3 8 ^ Bo 
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Using Lemma 4 . 2.5 and Lemma 4 . 2.6 
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objftct i’iitict iotj. C C'eepeetivfely t'eeofistir-iieted fttaetioti 3 value 
at i,j-th pixel, and MP is the total number o£ pixels used in 
summat X on . The first 50 images of sample 1 were used to 
estimate ot^’s. From (5.1.2) we see that the estimation 
problem reduces to that of finding a solution of the 
following system of linear equation 


e = A ^ (5.3.1) 

where 

-50x1 “ ^*1’ '*50^' ®k “ k-th image. 

-50x42 ~ ^50^’ * matrix consisting of 


covariance vectors as it’s rows and *42x1 vector of 
co-efficients to be estimated. Next to illustrate the case of 
unknown covariance kernel, the covariances were estimated 
from reconstructed images at the same 42 distances. These 
estimated covariance vectors Cj,’s substituted for C^'s in 
(5.3.1) give rise to the system of linear equations 

e = A (S (5.3.2) 

where Q is the 42x1 vector of co-efficients to be estimated. 

In order to solve the sysytem (5.3.1) [or (5.3.2) ] of 
linear equations the Residual Projection method was employed. 
This is because the columns of A are as such not independent, 
and so usual least square method to estimate the parameters 
is not applicable whereas the Residual Projection method 
gives a least squares solution. A brief description of the 
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TABLE 5.4 : SUMMARY OF ERRORS IK RECONSTRUCTION 
OF RANDOM IMAGES FROM SAMPLE 2 


PIC 

R 

c 

ERl 

BRS 

EIMaz 

E 

1 

16 

1.36925 

3.58375 

12 . 76621 

3 . 68718 


32 

0.49455 

0.39314 

6.19933 

0.59691 

2 

16 

1.37991 

3.78235 

16.52750 

9.55461 


32 

0.49257 

0.39152 

5.52099 

0.57670 

3 

16 

1.31931 

3.30063 

11.82067 

3.40341 


32 

0.47929 

0.36455 

6.39629 

0.55997 

4 

16 

1.33757 

3.39126 

11.37459 

3.48169 


32 

0.48486 

0.36446 

5.67.514 

0.56092 

5 

16 

1.32848 

3.41509 

10.49149 

3.47412 


32 

0.48738 

0.37286 

5.40592 

0.56061 

6 

16 

1.36158 

3.47953 

12.34249 

3.55429 


32 

0.48446 

0.37690 

6.15375 

0.57260 

7 

16 

1.34214 

3.40325 

11.06725 

3.49621 


32 

0.47869 

0.36821 

5.88351 

0-54924 

8 

16 

1.31688 

3.33913 

11.12923 

3.42109 


32 

0.46504 

0.34267 

5.62610 

0.51768 
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TABLE S. 12 : SUMMARY OF NORMALIZED ERRORS IN RECONSTRUCTION 
OF RANDOM IMAGES FROM SAMPLE 2 


PIC 

R 

c 

ERIN 

ERSN 

EN 

1 

16 

0.08940 

0.12034 

0.24074 


32 

0.03229 

0.04532 

0.03897 

2 

16 

0.08917 

0.12239 

0.61742 


32 

0.03183 

0.04477 

0.03727 

3 

16 

0.08390 

0.11284 

0.21644 


32 

0.03048 

0.04264 

0.03561 

4 

16 

0.08687 

0.11672 

0.22612 


32 

0.03149 

0.04351 

0.03643 

5 

16 

0.08578 

0.11629 

0.22432 


32 

0.03147 

0.04369 

0.03620 

6 

16 

0.08825 

'0.11756 

0.23037 


32 

0.03140 

0.04399 

0.03711 

7 

16 

0.08731 

0.11692 



32 

0.03114 

0.04373 

0.03573 

8 

16 

0.08736 

0.11819 

0-22695 


32 

0.03.085 

0.04305 

0.03434 
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Fig. 5 .2 Original siinulated random picture of sample 2 
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Then the estimate q of q is given bv 


q = Q R 

f' n n 


PROPOSITION: q 


q , w.p, 1 where q 

*pl opt 


<i R and 


the 


Moor e-Penrose g-inverse of Q. 


PROOF : Ue note that, at least along a subsequence of 
matrices, 

U >. U and V ». V , n »■ co .... 

where U and V are unitary. Since 


. in} 

A., 

j 


w.p. 1 

i = 1. . ,r 


1 0 

w.p. 1 

i = r+1 , . . 


wher e 

Xj.’s are singular 

values of 

<2, we have 


Q 

ri 

► Q* 

w.p. 1 

where 

the Moore-Penrose 

g-inverae 

Q* is 


1 





X 


1 

1 0 

I 

i 

t 

U . 



X 




0 


1 ® 



+ 

Q. is 


unitary 


( 6 . 3 . 2 ) 


This implies that along the subsequence 
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TABLE 6. S: ERROR IN CAT LINE INTEGRAL (NOISY) CASE 


PIC/FILTER OPT 

RAM 

RKS 

SHE 

CSI 

LIN 

4.199 

4.876 

5.399 

4.625 

4.509 

FAC 

3.847 

5.520 

5.999 

5.260 

5.194 

BRA 

5.052 

5.426 

5.995 

5.177 

5.035 

STA 

3.735 

4.371 

4.911 

4.095 

3.903 

COR 

13.924 

14.290 

15.000 

14.038 

14.038 

SAT 

4.320 

6.104 

6.606 

5.878 

5.853 

THO 

4.402 

4.938 

5.506 

4.662 

4.532 

SMO 

1.982 

3.094 

3.764 

2.694 

2.324 

PLA 

3.475 

4.247 

4.844 

3.931 

3.721 

BOY 

4.336 

4.492 

5.039 

4.238 

4.101 

MOM 

3.221 

4.382 

4.886 

4.119 

3.988 

UOR 

3 .639 

4.127 

4.657 

3.867 

3.670 

CHA 

5.896 

7.150 

7.711 

6.910 

6.884 

FIN 

3.327 

4.537 

5.021 

4.315 

4.215 

CRO 

4.083 

4.617 

5.133 

4.367 

4.195 

JET 

4.479 

5.023 

5.510 

4.760 

4.666 

GEO 

3.178 

4.230 

4.743 

3.985 

3.887 


AVERAGE 4.535 


5.378 5.925 5.113 4.983 




tion of Monalisa in CAT line integral noise- 
(a) CSI filter, (b) RAl'l filter, (c) Rx^S filte 
filter (e) original picture, (f) optimal 






stiruction of Chsi'sctsr's in < 
_f2rG6 CHSS (3) CSI filf^^/ 
KS filter, td) Shepp filter 
re (f) optimal filter- 














Reconstruction -or race oi. 
noisy cese CSI filter 
filter/ id) Shepp filter, 
(f) optimal filter. 
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Reconstruction- of Statue , of Liberty in PET sti 
integral noisy case (a) CS I f ilter,. (b) RAI“i fj 
(cl RKS'filtery- -Cdl filter, (e) original 

picture, (f) optimal; filter. 



